{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "#General import\n",
    "import numpy as np\n",
    "from scipy.integrate import ode\n",
    "from scipy.integrate import trapz\n",
    "import matplotlib.pyplot as plt\n",
    "import control\n",
    "import control.matlab\n",
    "#import seaborn as sns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "#sns.set()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD8CAYAAACMwORRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XtwnNWZ5/Hv07pa1s22ZEm+G7CNDTYGDDgXEmpIBocQqErIQnYrCRmIc5slmZ1JbZKtZXbyx9RkZiuZSUgl44FUyGUJCclSJAtDIKECTIJBGNtgG4wx4AvGlm+yJdnd6u5n/+huI8uSuyV199vvq9+nqou+HHU/r9v8dHze855j7o6IiERLLOgCRESk+BTuIiIRpHAXEYkghbuISAQp3EVEIkjhLiISQQp3EZEIUriLiESQwl1EJIKqg/rgtrY2X7BgQVAfLyISSs8999xBd2/P1y6wcF+wYAHd3d1BfbyISCiZ2RuFtNOwjIhIBCncRUQiSOEuIhJBCncRkQhSuIuIRJDCXUQkghTuIiIRFPlwTyTT/OTpN9hzZCDoUkREyiawi5jKoT+e5LZ7uvnTzkO0TKnhXz9+KavPmRF0WSIiJRfpnvv/Wb+LP+08xP+4dinTGmr42q9eIJXWhuAiEn2RDvcHNu7lojktfPo95/Dla85n58F+Hn5xX9BliYiUXGTDfceB42x58xjXr5wNwJoLOzm3fSrrntgZcGUiIqUX2XD/f5vfImbwoRVdAFTFjJsum8vmPb3sPqyTqyISbZEN9+d3H2FxRxMzm+tPPXfNBZ0APLLlraDKEhEpi7zhbmb1ZvaMmW0ysy1m9ncjtLnFzHrMbGP2dltpyi3cljePccGsltOemz9jKud3NincRSTyCum5x4E/c/eLgJXAGjNbPUK7+9x9ZfZ2V1GrHKMDx0/SczzOBbOaz3jtmgs66X7jCIf7EwFUJiJSHnnD3TP6sg9rsreKnk+45c1jACOG+3sWt+EOz7x2qNxliYiUTUFj7mZWZWYbgQPAo+6+foRmHzGzzWZ2v5nNLWqVY7Q1G+7LRgj35bNbmVJTxdM7D5e7LBGRsiko3N095e4rgTnA5WZ24bAmvwYWuPsK4FHgnpHex8zWmlm3mXX39PRMpO6z2vJmL/NnNNBUX3PGa7XVMVYtmMafXlXPXUSia0yzZdz9KPA4sGbY84fcPZ59eBdw6Sg/v87dV7n7qvb2vPu7jtuOA30s7mga9fXV58zg5f3HOdQXH7WNiEiYFTJbpt3MWrP3pwDvB14a1qZryMPrgW3FLHIs3J09R04wd1rDqG1WnzMdgGdfP1KuskREyqqQnnsX8LiZbQaeJTPm/hsz+7qZXZ9tc3t2muQm4HbgltKUm9/RgUEGEinmTJsyapsLZrVQU2Vs3H20jJWJiJRP3lUh3X0zcPEIz98x5P5Xga8Wt7Tx2XPkBACzzxLu9TVVLOtq5vld6rmLSDRF7grV3LrtZ+u5A1w8bxov7O3VKpEiEkkRDPdMz33OWcbcAVbObWUgkWL7/uPlKEtEpKwiGO4DNNVV0zLlzGmQQ62c2wrA87s07i4i0RO5cN979MRZx9tz5s9ooGVKDS++2VuGqkREyity4b7nyIm8QzIAZsayruZTSxWIiERJpMI9N8c938nUnGWzmnlp3zGSqXSJKxMRKa9IhftAIkVfPElnS33+xmQWFosn07x+qL/ElYmIlFekwv1gdjmBtsa6gtrnFhbT0IyIRE3Ewj2zRvuMxtqC2p/b3khtVezUKpIiIlERqXDPLQTWNrWwnntNVYzFnY1s3adwF5FoiVS453rubU2F9dwBLuhqYeubx3DXlaoiEh2RCvdcz3361MLDfdmsZg71J9h/TMv/ikh0RCrcD/bFaa6vpq66quCfyZ1U3bpPFzOJSHREK9z7EwXPlMk5vzOzqYdOqopIlEQq3A/1xccc7k31NSyY0aDpkCISKZEK94N9iYKnQQ61bJaWIRCRaIlUuB/qi48r3Jd2NrPr8AADiWQJqhIRKb/IhPtgKs2RgcExD8sALMpupv3K/r5ilyUiEohCNsiuN7NnzGxTdp/UvxuhTZ2Z3WdmO8xsvZktKEWxZ3OkP3d16tjDfUn2pOrL2rhDRCKikJ57HPgzd78IWAmsMbPVw9rcChxx9/OAbwHfKG6Z+fVk57i3j2NYZt70BuqqY2x/S+EuItGQN9w9IzdeUZO9Db+c8wbgnuz9+4GrzcyKVmUBjg4MAtDaMPZwr4oZizoa1XMXkcgoaMzdzKrMbCNwAHjU3dcPazIb2A3g7kmgF5hRzELzOXYiE+75ttcbzeKZTRpzF5HIKCjc3T3l7iuBOcDlZnbheD7MzNaaWbeZdff09IznLUbVmw335vGGe2cTbx07SW/2XwAiImE2ptky7n4UeBxYM+ylvcBcADOrBlqAQyP8/Dp3X+Xuq9rb28dX8SiOnZxYz31JdsbM9gMamhGR8Ctktky7mbVm708B3g+8NKzZg8Ans/dvBH7vZV5msffEIFUxY2pt4evKDLU4N2NGJ1VFJAKqC2jTBdxjZlVkfhn83N1/Y2ZfB7rd/UHgbuDHZrYDOAzcXLKKR3HsRJLm+mrGex53Vks9jXXVbNdJVRGJgLzh7u6bgYtHeP6OIfdPAh8tbmlj03ticNxDMgBmmRkzCncRiYLIXKHae2Jw3CdTc5Z0NPHyW8e1cYeIhF5kwv3YyYn13AEWdzRxZGDw1I5OIiJhFZlw7z0xSHP9BHvu2ZOqGpoRkbCLTLgfO5Gc8LDMoo5GQDNmRCT8IhHu7s6xCZ5QBWhvrGNaQw2vaK67iIRcJMI9nkyTSKVpnlLIzM7RmRmLsydVRUTCLBLh3jvBdWWGWtLZxPb9fZoxIyKhFqlwn+gJVchs3NEXT/Jm78kJv5eISFAiEe4TXRFyqFNrzGhoRkRCLBLhXsxhmcXZGTOaDikiYRaJcM+tCDnRqZCQ2eyjo7lOJ1VFJNQiEe65Ndib6yc2WyZncUeTlv4VkVCLRLj3xZMANBXhhCpkxt1f2d9HKq0ZMyISTpEI9/5EitqqGLXVxTmcxZ1NxJNpdh0eKMr7iYiUWzTCPZ6koW58m3SMJDdjRuPuIhJWkQj3vniSqbXFGW+Ht9eY0YwZEQmrSIT7QDxFY13xwr2htpp50xt4WeEuIiEViXDvTxR3WAayM2Y0LCMiIVXIBtlzzexxM9tqZlvM7IsjtLnKzHrNbGP2dsdI71UqffFkUXvuAEs6G3ntYD/xZKqo7ysiUg6FJGIS+Gt332BmTcBzZvaou28d1u5Jd7+u+CXm1x9P0tFUX9T3XNzRRDLtvHawn/M7m4v63iIipZa35+7u+9x9Q/b+cWAbMLvUhY1FfzxV9GGZ3K5MmjEjImE0pjF3M1sAXAysH+Hld5jZJjN72MwuKEJtBetPFH9Y5py2RqpjphkzIhJKBSeimTUCvwS+5O7Hhr28AZjv7n1mdi3wALBohPdYC6wFmDdv3riLHq4/nmRqkcO9tjrGwrapvPxWX1HfV0SkHArquZtZDZlg/6m7/2r46+5+zN37svcfAmrMrG2EduvcfZW7r2pvb59g6RmJZJrBlDO1trjDMpC5UlU9dxEJo0JmyxhwN7DN3b85SpvObDvM7PLs+x4qZqGj6c+uK1PsnjtkrlTddXiAgUSy6O8tIlJKhSTiu4CPAy+Y2cbsc18D5gG4+/eBG4HPmVkSOAHc7GXap66vhOG+OLsMwSv7+7hobmvR319EpFTyJqK7PwVYnjZ3AncWq6ixGEhk5qEXc/mBnFMzZvYfV7iLSKiE/grVt3vuxR9znze9gbrqmK5UFZHQCX24l3LMvSpmLOpo1BozIhI6oQ/33MnOUgzLQHaNGYW7iIRM6MO9L54Zcy/2RUw5Szqa2H8sztGBREneX0SkFEIf7rlhmWIvP5CzOHtSdft+XcwkIuER+nDPnVAtZc8d0Li7iIRK6MN9IJGkKmbUFWn/1OG6WuppqqvWjBkRCZXQh3t/PEVDbRXZC2SLzsxY3NmknruIhErow70UG3UMl5sxU6aLbkVEJiz04T6QSNJQgkXDhlrS0cjRgUF6+uIl/RwRkWIJfbifSKRoKNEc95xTM2a0/K+IhET4w30wxZSaUvfcM+H+0lvDl7EXEalMEQj3NPUlHpaZ0VhHW2MdL2nGjIiEROjD/WQixZSa0h/G0q4mtu1Tz11EwiH04V6OYRmAZbOaeWV/H4OpdMk/S0RkoqIR7iUelgFY1tVMIpXm1R6dVBWRyhf6cD+ZSFFfhp770q5mALa+qaEZEal8oQ/3cg3LnNM2ldrqmMbdRSQUCtkge66ZPW5mW81si5l9cYQ2ZmbfNrMdZrbZzC4pTbmnG0ylSaa9LOFeXRVjSUcT2/ZpxoyIVL5Ceu5J4K/dfRmwGviCmS0b1uYDwKLsbS3wvaJWOYqTg5m13Msx5g6ZGTNb9x3TMgQiUvHyhru773P3Ddn7x4FtwOxhzW4AfuQZTwOtZtZV9GqHOZEN93KMuUNm3P1wf4IDx7UMgYhUtjGNuZvZAuBiYP2wl2YDu4c83sOZvwAws7Vm1m1m3T09PWOrdAQnE5lpieUYloHMjBmArRp3F5EKV3C4m1kj8EvgS+4+rnRz93XuvsrdV7W3t4/nLU5zoszDMudrxoyIhERB4W5mNWSC/afu/qsRmuwF5g55PCf7XEmdCvcy9dxbptQwu3WKZsyISMUrZLaMAXcD29z9m6M0exD4RHbWzGqg1933FbHOEZ1IlHfMHTJXqmpYRkQqXSFr5b4L+DjwgpltzD73NWAegLt/H3gIuBbYAQwAnyp+qWcq92wZyIy7P7ZtP/3xJFNLvEmIiMh45U0nd38KOOsedp6ZG/iFYhVVqLdny5TvWqzls1twz5xUvWzB9LJ9rojIWIT6CtXcsEy5xtwBls9pAeCFPb1l+0wRkbEKd7iX+YQqQEdzPe1Ndby4V+EuIpUr1OGeG3Mv9WYdwy2f3cILCncRqWChDvcghmUALpzdwqs9fQwkkmX9XBGRQoU73AdTVMeMmqryHsby2S2kXRcziUjlCn24l7vXDplwBzQ0IyIVK9ThfnIwVfbxdoCO5syG2Qp3EalUoQ73E4lgeu5mxvLZzZoxIyIVK9zhHtCwDGSGZnYc0ElVEalMIQ/3dCDDMpCZMZN2tIiYiFSkUIf7yUSKKWVcemCoFXNaAdisK1VFpAKFOtyDHJbpbKmns7mejbuPBvL5IiJnE/5wD2hYBmDl3Fae36VwF5HKE+pwPzmYor46uHC/eF4ruw4PcKhPe6qKSGUJebgHd0IV4OJ50wA0NCMiFSfU4R5PpqirDu4Qls9uoSpmGpoRkYoT8nBPUxfgsMyU2irO72zi+d1HAqtBRGQkoQ13dyeRTAfac4fMSdVNu3tJpT3QOkREhipkg+wfmNkBM3txlNevMrNeM9uYvd1R/DLPFE+mAagLaJ57zsXzptEXT/JqT1+gdYiIDFVIMv4QWJOnzZPuvjJ7+/rEy8rvVLgHOCwDmRkzAM/v0tCMiFSOvOHu7k8Ah8tQy5jEk5mNOoIellk4YyrN9dWaMSMiFaVYyfgOM9tkZg+b2QWjNTKztWbWbWbdPT09E/rA+GCm514bcLjHYsbKedPY8IbCXUQqRzGScQMw390vAr4DPDBaQ3df5+6r3H1Ve3v7hD707WGZ4M8JXzZ/Gi/vP87RgUTQpYiIAEUId3c/5u592fsPATVm1jbhyvJ4e1gm2DF3gMsWTgeg+3WNu4tIZZhwuJtZp5lZ9v7l2fc8NNH3zadSZstAZjpkbVWMZ16vuFMTIjJJVedrYGb3AlcBbWa2B/hboAbA3b8P3Ah8zsySwAngZncv+aTv3Jh7JQzL1NdUsWJOC8+8pnAXkcqQN9zd/WN5Xr8TuLNoFRWokoZlAC5fOJ11T+ykP55kal3eP1YRkZIKvts7TpV0QhUy4Z5Mu9aZEZGKUBnJOA65cK+vgDF3gEvnTyNmaNxdRCpCZSTjOCQq5ArVnKb6GpbNauaZ10p+LllEJK/QhnulXKE61GULpvP8rqOnfvGIiASlcpJxjN6eLVMZPXeAKxZOJ55Ms3mPxt1FJFjhDfcKmueec8XCGZjBH1/V0IyIBKtyknGMcsMytVWVcwjTptZywaxmnnrlYNCliMgkVznJOEbxZJraqhixmAVdymnefV47G3YdoS+eDLoUEZnEwhvug8HvwjSSKxe1kUy7Zs2ISKAqLx0LFE+mKmq8PefS+dOoq47xpIZmRCRAlZeOBQp6c+zR1NdUcfnC6fzHDoW7iAQn5OFemeW/+7w2tu/vY/+xk0GXIiKTVGWmYwHig6nAd2EazbsXZZaz16wZEQlKZaZjAeLJNHU1lTcsA7C0s5kZU2t5SkMzIhKQEId7qmKHZWIx48pFbfxhew+pdMmXthcROUNlpmMBKnnMHeDqpR0c7k/w/C5tvSci5Ve56ZhHpc5zz3nvknaqY8Zj2w4EXYqITEKVm455ZIZlKnPMHaC5vobLF07nd9v2B12KiExCecPdzH5gZgfM7MVRXjcz+7aZ7TCzzWZ2SfHLPFOlD8tAZmjmlQN97Do0EHQpIjLJFJKOPwTWnOX1DwCLsre1wPcmXlZ+mdkylR3u71s6E4DH1HsXkTLLm47u/gRwtr3jbgB+5BlPA61m1lWsAkcTH6zsYRmA+TOmct7MRn73ksJdRMqrGF3f2cDuIY/3ZJ87g5mtNbNuM+vu6emZ0IeGYVgG4OqlM1m/8zDHTg4GXYqITCJlTUd3X+fuq9x9VXt7+0Teh0QqHOH+58s6Saadx7aq9y4i5VOMdNwLzB3yeE72uZIZTDnuVOwVqkNdMq+V2a1T+M3mfUGXIiKTSDHC/UHgE9lZM6uBXncvaZJV4ubYozEzPriiiydf6eHoQCLockRkkihkKuS9wJ+AJWa2x8xuNbPPmtlns00eAnYCO4B/Az5fsmqzTu2fGoJwB7huRReDKeeRLW8FXYqITBLV+Rq4+8fyvO7AF4pWUQFy4V6pq0IOt3x2C/NnNPCbzfu46bJ5QZcjIpNAONJxmETIwt3MuG5FF3989RCH+uJBlyMik0A40nGYwVQm3GuqwlP+dStmkUo7D7+ooRkRKb3wpOMQp3ruIQr38zubWDSzkf/7fEknEomIAGEN91zPPSTDMpAZmvnoqjk898YRdhw4HnQ5IhJx4UnHIQZzs2VC1HMH+PAlc6iOGfc9uzt/YxGRCQhXOmaFsecO0NZYx/uWdvCrDXtPDS2JiJRCuNIxK3dCNUxj7jk3XTaXQ/0Jfq/FxESkhMKXjrx9QjVMs2Vy3rO4nc7meg3NiEhJhS8dgUQqs+l0bbUFXMnYVcWMGy+dwx+297D7sDbxEJHSCGe4n5oKWfkLh43kP18xDzPjnj++HnQpIhJRoQz3UxcxhbDnDjCrdQofXN7Ffc/u5rjWeReREgh1uIfxhGrOre9eyPF4kp937wm6FBGJoFCm46kTqiGbCjnURXNbWTV/Gj/842uk0h50OSISMaFMx0QEeu4At125kN2HT/BbLQUsIkUWynQM49oyI3n/sk4WzGjgO7/fQWblZBGR4ghlOg6m0lTHjFgsnCdUc6pixu1XL2LrvmM8skUXNYlI8YQy3BPJdCgvYBrJ9RfN4pz2qfzzY9tJa+xdRIqkoIQ0szVm9rKZ7TCzr4zw+i1m1mNmG7O324pf6tsGU05NVbh77TnVVTG+ePUiXnrruNZ6F5GiKWQP1Srgu8AHgGXAx8xs2QhN73P3ldnbXUWu8zTxZJra6nBewDSS61bMYtHMRr712HaSKS0oJiITV0jP/XJgh7vvdPcE8DPghtKWdXaDqTS1Eem5Q2bs/W+uWcKOA3385Ok3gi5HRCKgkHCfDQxd5WpP9rnhPmJmm83sfjObW5TqRjGYSodm/9RC/fmyDq5c1MY3H93O4f5E0OWISMgVKyF/DSxw9xXAo8A9IzUys7Vm1m1m3T09PeP+sCidUM0xM+64bhn9iRT/+7cvB12OiIRcIQm5FxjaE5+Tfe4Udz/k7vHsw7uAS0d6I3df5+6r3H1Ve3v7eOoFotlzB1jU0cQn37GAe5/ZxQt7eoMuR0RCrJCEfBZYZGYLzawWuBl4cGgDM+sa8vB6YFvxSjxTPII995wvvm8RbY11fPn+TdqtSUTGLW9CunsS+EvgETKh/XN332JmXzez67PNbjezLWa2CbgduKVUBUPuhGo0w71lSg3/8OHlvPTWcb79u1eCLkdEQqq6kEbu/hDw0LDn7hhy/6vAV4tb2ugSyTQNtQWVHkpXL+3gxkvn8L0/vMr7l3Vw0dzWoEsSkZAJZfc3ShcxjeaODy1jZlMdf3XfRq35LiJjFspwTySjeUJ1qOb6Gr5100reODzAl3+xWQuLiciYhDIhB1PRPaE61OpzZvCVNefz71ve4l+f2Bl0OSISIqFMyEREp0KO5LYrF/LB5V3847+/xOMvHQi6HBEJiVAmZCIZ3dkyw5kZ37hxBUu7mvn8TzewYdeRoEsSkRAIZUJG9SKm0TTWVfPDT13OzOY6/uKHz7LjwPGgSxKRChfKhIzi8gP5tDfV8eO/uILqWIyP/dt6Xn5LAS8iowtlQmamQoay9AmZN6OBez99BTGDm9b9iU27jwZdkohUqNAlpLtPqhOqwy3qaOIXn3knTfXV/Je71vP7l7Q9n4icKXQJOZjKzPeO0nruYzVvRgO/+Mw7WdDWwK33dPPdx7XBtoicLoThnllMa7L23HM6W+r5xWfeyYdWzOKfHnmZT//oOQ72xfP/oIhMCqFLyNxKiZNxzH24KbVV/MvNK/mf1y3jie09XPOtJ/jtFu3DKiIhDHf13E9nZtz67oX8+r++m47metb++Dluu6ebXYcGgi5NRAIUuoSMq+c+oiWdTTzwhXfxlQ+czx9fPcj7vvUH/v6hbRqqEZmkQpeQuZ57nXruZ6itjvHZ957L439zFdet6OKuJ3dy5Tce5+8f2saeI+rJi0wmoUvIREo993w6muv55n9ayaP/7b1cc0EHdz/1Gu/5x8f5zI+7+d22/drhSWQSCN2OF4PJzJQ/hXt+57Y38s83X8yX15zPT55+g589s4tHtuyntaGGD1zYxYcu6mLV/Ok6fyESQaEL90QqBeiE6ljMbp3Cf19zPn/1vsU8+UoPD256kwee38u9z+xiam0Vq8+ZwZWL2lh97gwWzWyiKjZ5ryEQiYqCwt3M1gD/AlQBd7n7Pwx7vQ74EXApcAi4yd1fL26pGYlTPXcF0FjVVse4emkHVy/tYCCR5IntB3lqRw9PvnKQ32WXE26orWL57BZWzm1lSWcT581s5Nz2RqbWha4fIDKp5f0/1syqgO8C7wf2AM+a2YPuvnVIs1uBI+5+npndDHwDuKkUBeuEanE01Faz5sJO1lzYCcCuQwN0v3GYTbuPsnH3UX7wH6+duhoYYFZLPefObGR26xS6WqbQ1VJPV2s9XS31TJ9aR8uUGvX4RSpIId2xy4Ed7r4TwMx+BtwADA33G4D/lb1/P3CnmZmX4Jp4XcRUGvNmNDBvRgMfvmQOkPlz3nW4nx0H+ni1J/PfnT19bNt3fMTplWbQOqWGaQ21TJtay7SGGqbWVdNQW83U2ioa6rL/ra3KPFdXRX1NFTVVsezNqKmKUVudeVwds1P3a6qM6lgMM6iKGTEzYpaZ4y8iIysk3GcDu4c83gNcMVobd0+aWS8wAzhYjCKH0kVM5VFbHeO8mU2cN7PpjNfiyRT7e+O82XuC/cdOcrg/wZH+BEcGBjk8kODoQIK9R08ykEgykEgxEE/Sn0iVpM5M2GeCPmYQM6PKDDOInfpFkGsDxtu/EHK/G4b+ihjtF8bQp0+7P8L7jfaep71znrYSbTdfNpfbrjynpJ9R1oFUM1sLrAWYN2/euN5jZnMd1y7vpGVKTTFLkzGoq6461dMvVDrtnEym6I+nGEgk6Y+nODGYIplKM5hyBlPp7C1zP5F9nDz1mpN2x91JO6TdSaeH3B/hOXcnlX3N3UkPmQHqZP5ROfTflkP/mXn68yM3Or29j/J84W3R2m+TRltjXck/o5Bw3wvMHfJ4Tva5kdrsMbNqoIXMidXTuPs6YB3AqlWrxvVX+dL507l0/vTx/KgEKBYzGmozwzRQ+r/YIpNdIWMbzwKLzGyhmdUCNwMPDmvzIPDJ7P0bgd+XYrxdREQKk7fnnh1D/0vgETJTIX/g7lvM7OtAt7s/CNwN/NjMdgCHyfwCEBGRgBQ05u7uDwEPDXvujiH3TwIfLW5pIiIyXppyIiISQQp3EZEIUriLiESQwl1EJIIU7iIiEWRBTUc3sx7gjXH+eBslWNogBCbjceuYJwcdc+Hmu3t7vkaBhftEmFm3u68Kuo5ym4zHrWOeHHTMxadhGRGRCFK4i4hEUFjDfV3QBQRkMh63jnly0DEXWSjH3EVE5OzC2nMXEZGzqOhwN7M1Zvayme0ws6+M8Hqdmd2XfX29mS0of5XFVcAx32JmPWa2MXu7LYg6i8nMfmBmB8zsxVFeNzP7dvbPZLOZXVLuGoutgGO+ysx6h3zPd4zULkzMbK6ZPW5mW81si5l9cYQ2kfquCzzm0nzXnt2xptJuZJYXfhU4B6gFNgHLhrX5PPD97P2bgfuCrrsMx3wLcGfQtRb5uN8DXAK8OMrr1wIPk9mRbjWwPuiay3DMVwG/CbrOIh9zF3BJ9n4TsH2Ev9+R+q4LPOaSfNeV3HM/tTG3uyeA3MbcQ90A3JO9fz9wtYV7E8pCjjly3P0JMvsAjOYG4Eee8TTQamZd5amuNAo45shx933uviF7/ziwjcz+y0NF6rsu8JhLopLDfaSNuYf/oZy2MTeQ25g7rAo5ZoCPZP/Jer+ZzR3h9agp9M8lat5hZpvM7GEzuyDoYoopO4R6MbB+2EuR/a7PcsxQgu+6ksNdRvZrYIG7rwAe5e1/uUi0bCBzmflFwHeABwKup2jMrBH4JfAldz8WdD3lkOeYS/JdV3K4j2Vjbs62MXeI5D1mdz/k7vHsw7uAS8tUW5AK+bsQKe5+zN37svcfAmrMrC3gsibMzGrIhNxP3f1tGtbdAAAA7ElEQVRXIzSJ3Hed75hL9V1XcrhPxo258x7zsPHH68mM4UXdg8AnsjMpVgO97r4v6KJKycw6c+ePzOxyMv+vhrnjQvZ47ga2ufs3R2kWqe+6kGMu1Xdd0B6qQfBJuDF3gcd8u5ldDyTJHPMtgRVcJGZ2L5kZA21mtgf4W6AGwN2/T2b/3muBHcAA8KlgKi2eAo75RuBzZpYETgA3h7zjAvAu4OPAC2a2Mfvc14B5ENnvupBjLsl3rStURUQiqJKHZUREZJwU7iIiEaRwFxGJIIW7iEgEKdxFRCJI4S4iEkEKdxGRCFK4i4hE0P8HDJM6UqB5YuwAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Define a continuos-time power spectrum with bandwidth 10, and \n",
    "# standard deviation 3 (when filtered by white noise)\n",
    "\n",
    "wu = 10 # bandwidth of the force disturbance\n",
    "std_du = 3\n",
    "Hu = control.TransferFunction([1], [1/wu, 1])\n",
    "Hu = Hu * Hu\n",
    "t_imp = np.arange(0,wu/4, 0.1/wu)\n",
    "t, y = control.impulse_response(Hu,t_imp)\n",
    "plt.plot(t,y)\n",
    "std_tmp = np.sqrt(trapz(y**2,t))\n",
    "Hu = Hu/(std_tmp)*std_du"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAD8CAYAAABpcuN4AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3X10XPV95/H3V6MHW7YsP8mykG1ksGwQT4YIA2loE5wcTJrGJCWJ2bMpaUlJdqFtNj09gW6XTeh6W3JOwiYpJIcWEkrSGJd0G8GScGggDaGJbUGMwTbGwgb8hC1btmVbD6MZffePe2WGYSSNpXmez+ucObpz7+/+5vfT2L+vfg/3XnN3RERERlOR7wKIiEhhU6AQEZExKVCIiMiYFChERGRMChQiIjImBQoRERmTAoWIiIxJgUJERMakQCEiImOqzHcBMmHu3Lne0tKS72KIiBSV559//rC7N4yXriQCRUtLC52dnfkuhohIUTGzN9JJp6EnEREZkwKFiIiMSYFCRETGpEAhIiJjUqAQEZExKVCIiMiYFChERGRMChRnYOv+43z756+xp6cv30UREcmZkrjgLhdeeauXj933H0Rjw/z9s7t4/E/ex1kzp+a7WCIiWaceRZrW/r/tTK+p5Id/fCX90Th//fi2fBdJRCQnFCjSsKenj2d3HuYP39vCVefO4bNXL+anW9+i69DJfBdNRCTrFCjS8OTWtwC4/tJmAG56bwuVFca6jW/ms1giIjmhQJGG57oOc87caSycXQvA3Ok1/M7SBh7fcoDhYc9z6UREskuBYhxD8WE27u7hqnPnvGP/711yFm/1DtD5xtE8lUxEJDfSChRmtsrMdphZl5ndnuJ4jZk9Eh7fYGYtCcfuCPfvMLNrw31TzGyjmb1oZlvN7CsJ6b9nZrvNbHP4Wj75ak5c16GTnIrGubxl9jv2f+C8eVRWGD/fcShPJRMRyY1xA4WZRYB7geuANuBGM2tLSnYzcNTdlwD3AHeH57YBa4ALgFXAfWF+g8A17n4JsBxYZWZXJuT3F+6+PHxtnlQNJ+mVt3oBaDtrxjv2z5hSxWVnz+LfX+3OR7FERHImnR7FCqDL3Xe5exRYB6xOSrMaeCjcfhRYaWYW7l/n7oPuvhvoAlZ4YGTJUFX4KsjB/u0HTlBdWcE5c6e969jvLG1g6/5eDp0YyEPJRERyI51A0QzsSXi/N9yXMo27x4DjwJyxzjWziJltBg4BT7n7hoR0a81si5ndY2Y1Z1CfjNt+oJeljdOpjLz7V3V161wAfr2rJ9fFEhHJmbxNZrt73N2XAwuAFWZ2YXjoDuA84HJgNvClVOeb2S1m1mlmnd3d2Rv+eeWtE5w3f0bKY21NM6itjtD5ugKFiJSudALFPmBhwvsF4b6UacysEqgHjqRzrrsfA54hmMPA3Q+EQ1ODwHcJhr7exd3vd/d2d29vaBj32eATcnIwRveJQc5pePewE0BlpILLFs1i0+ta+SQipSudQLEJaDWzxWZWTTA53ZGUpgO4Kdy+AXja3T3cvyZcFbUYaAU2mlmDmc0EMLOpwIeAV8L3TeFPA64HXp5MBSdj5OZ/Z89OHSgALm+ZzStv9dI7MJSrYomI5NS4NwV095iZ3QY8CUSAB919q5ndBXS6ewfwAPCwmXUBPQTBhDDdemAbEANudfd4GAweCldAVQDr3f3x8CN/YGYNgAGbgc9nssJn4s0wUCwKL7RL5fKWWbjD828c5QPL5uWqaCIiOZPW3WPd/QngiaR9dyZsDwCfGOXctcDapH1bgEtHSX9NOmXKhTePjB8oli+aSYXBCwoUIlKidGX2GN7s6aN+ahX1tVWjpqmtrmRpYx0v7Tuew5KJiOSOAsUY3uzpG7M3MeKi5npe2nucYFpGRKS0KFCMYf+xfprTeDjRxQvqOXIqyr5j/TkolYhIbilQjOFg7wDz66eMm+7iBTMBeGmvhp9EpPQoUIyiPxqndyBGQ934F4af11RHVcTYonkKESlBChSjGLl/U+OM8XsUNZURls2vU49CREqSAsUoDvYOAtA4I71bTV3UPJMte49pQltESo4CxSgO9qbfowC4sHkGvQMx9h7VhLaIlBYFilGcDhR16QWK85uCGwduP9CbtTKJiOSDAsUoDp0YpKayghlT07p4nWWNdZgFz68QESklChSjONg7QOOMKQT3JhzftJpKzp5de/qJeCIipUKBYhQHeweYl8bS2ETnN83Q0JOIlBwFilH0nIoyd/qZB4o3evo4NRjLUqlERHJPgWIUPaeGmDWt+ozOOW9+He6w46DmKUSkdChQpODuHO2LMnva6HeNTWVk5dMrmtAWkRKiQJFC70CM+LAzq/bMehQLZk2lrqZS8xQiUlIUKFI4eioKwOwzHHoyM85rqlOgEJGSokCRQk9fECjOdI4C4Lz5M3jlrRO6lYeIlIy0AoWZrTKzHWbWZWa3pzheY2aPhMc3mFlLwrE7wv07zOzacN8UM9toZi+a2VYz+0pC+sVhHl1hnmfeWk/S6R7FGQ49QTBPcXJQt/IQkdIxbqAwswhwL3Ad0AbcaGZtScluBo66+xLgHuDu8Nw2YA1wAbAKuC/MbxC4xt0vAZYDq8zsyjCvu4F7wryOhnnnVM8Eh54Als2fDsDOQ5rQFpHSkE6PYgXQ5e673D0KrANWJ6VZDTwUbj8KrLTgkubVwDp3H3T33UAXsMIDJ8P0VeHLw3OuCfMgzPP6CdZtwo5OYuhpybw6AHa8dXKclCIixSGdQNEM7El4vzfclzKNu8eA48Ccsc41s4iZbQYOAU+5+4bwnGNhHqN9Vtb1nBqiOlLBtOrIGZ9bP7WK+TOmsFPXUohIicjbZLa7x919ObAAWGFmF57J+WZ2i5l1mllnd3d3Rst29FSUWdOq0r7PU7Kl8+t00Z2IlIx0AsU+YGHC+wXhvpRpzKwSqAeOpHOuux8DniGYwzgCzAzzGO2zRs67393b3b29oaEhjWqkr6cvesbXUCRaOm86XYdOEh/WyicRKX7pBIpNQGu4GqmaYHK6IylNB3BTuH0D8LQH60M7gDXhqqjFQCuw0cwazGwmgJlNBT4EvBKe80yYB2GeP5549SbmeN8QM2vP7KrsREvn1zEYG+bNnr4MlkpEJD/GDRThfMFtwJPAdmC9u281s7vM7KNhsgeAOWbWBXwRuD08dyuwHtgG/BS41d3jQBPwjJltIQhET7n742FeXwK+GOY1J8w7p3oHhpgxZRKBojGY0H5Vw08iUgLSeiqPuz8BPJG0786E7QHgE6OcuxZYm7RvC3DpKOl3Eay0ypve/iFmTJ14oGidFy6RPXiCay+Yn6liiYjkha7MTuHEQIy6Kek92S6VaTWVLJg1lR0HtURWRIqfAkWS+LBzYjA2qaEnCIaftERWREqBAkWSk+FDhybTowBobZzOa90nGYoPZ6JYIiJ5o0CRpLd/CGBScxQAyxrrGIo7bxw5lYliiYjkjQJFkt6BMFBkYOgJ4FXNU4hIkVOgSHJiIBh6mjHJoacl86ZjBjve0jyFiBQ3BYokmRp6mlIV4ezZtbqLrIgUPQWKJL0DmZnMBmhtrFOPQkSKngJFkhMZmqOAYEL79SN9DMbik85LRCRfFCiS9PZnskcxnfiws/uwVj6JSPFSoEjSOzBEbXWEysjkfzUjK580/CQixUyBIsmJSd4QMNE5DdOIVBg7tURWRIqYAkWS3v7J3ecpUU1lhJY5tXqIkYgUNQWKJL0Dk7tzbLJl83XPJxEpbgoUSU4MxCZ9sV2ipY11vNHTR39UK59EpDgpUCQ5ORhjWk3mAsWyxjrcoeuQ5ilEpDgpUCQ5ORhjegYDxdL54conDT+JSJFSoEjSl+Eexdmza6murNBjUUWkaKUVKMxslZntMLMuM7s9xfEaM3skPL7BzFoSjt0R7t9hZteG+xaa2TNmts3MtprZnyWk/7KZ7TOzzeHrw5OvZnqGh52+oTjTqiMZy7MyUsGShum6lkJEita4gcLMIsC9wHVAG3CjmbUlJbsZOOruS4B7gLvDc9uANcAFwCrgvjC/GPDn7t4GXAncmpTnPe6+PHy941nd2dQ/FMedjPYoIFj5pB6FiBSrdHoUK4Aud9/l7lFgHbA6Kc1q4KFw+1FgpZlZuH+duw+6+26gC1jh7gfc/QUAdz8BbAeaJ1+dyTkVDW7fUZvhQLG0sY4Dxwc4Ht6ZVkSkmKQTKJqBPQnv9/LuRv10GnePAceBOemcGw5TXQpsSNh9m5ltMbMHzWxWGmXMiL7BYAlrJoeeAJbNnw6g6ylEpCjldTLbzKYDPwK+4O694e5vA+cCy4EDwNdGOfcWM+s0s87u7u6MlGfkedmZHno6fc8nBQoRKULpBIp9wMKE9wvCfSnTmFklUA8cGetcM6siCBI/cPd/GUng7gfdPe7uw8DfEwx9vYu73+/u7e7e3tDQkEY1xtcXHelRZDZQNM+cyrTqCK9qQltEilA6gWIT0Gpmi82smmByuiMpTQdwU7h9A/C0u3u4f024Kmox0ApsDOcvHgC2u/vXEzMys6aEtx8DXj7TSk3U23MUmR16MjOWzq9Tj0JEitK4fzq7e8zMbgOeBCLAg+6+1czuAjrdvYOg0X/YzLqAHoJgQphuPbCNYKXTre4eN7P3AZ8GXjKzzeFH/WW4wumrZrYccOB14HMZrO+YRuYoMnnB3YhljXU8te1gxvMVEcm2tFrEsAF/ImnfnQnbA8AnRjl3LbA2ad8vARsl/afTKVM2nArnKGozPJkNwTzFuk17OHxykLnTazKev4hItujK7AQjQ0+ZnqOAtye0NU8hIsVGgSLB6cnsLAw9LQ2XyGqeQkSKjQJFgpODMaoiRnVl5n8tDdNrmFVbpSu0RaToKFAk6BuMUZuFYScIVz411umeTyJSdBQoEpyKxrOy4mlEcM+nkwQrh0VEioMCRYJTg7GsrHgasbSxjpODMfYfH8jaZ4iIZJoCRYJT0XhWJrJHLJuvlU8iUnwUKBIEDy3KYo9inu75JCLFR4EiwcksTmYD1NdWMX/GFPUoRKSoKFAk6MvyZDagez6JSNFRoEjQF40xNYuT2QDLGqez89BJ4sNa+SQixUGBIkFfNLPPy05laWMd0dgwbxw5ldXPERHJFAWKkLvTPxRnalWWexThyqdXNE8hIkVCgSI0GBvGHabkoEcRqTC2H+gdP7GISAFQoAgNDAU3BMx2j2JKVYRz5k5ToBCRoqFAEerPUaAAaDtrBtv2K1CISHFQoAgNDA0DZH3VE8D5TTPYf3yAY33RrH+WiMhkKVCE+sNnUdRU5qBH0TQDgG0afhKRIpBWoDCzVWa2w8y6zOz2FMdrzOyR8PgGM2tJOHZHuH+HmV0b7ltoZs+Y2TYz22pmf5aQfraZPWVmO8OfsyZfzfGdHnrKUY8C0PCTiBSFcQOFmUWAe4HrgDbgRjNrS0p2M3DU3ZcA9wB3h+e2AWuAC4BVwH1hfjHgz929DbgSuDUhz9uBn7l7K/Cz8H3W5WoyG6ChroaGuhq2H9ASWREpfOn0KFYAXe6+y92jwDpgdVKa1cBD4fajwEozs3D/OncfdPfdQBewwt0PuPsLAO5+AtgONKfI6yHg+olV7cyMDD3lIlBAMPykoScRKQbpBIpmYE/C+7283ai/K427x4DjwJx0zg2HqS4FNoS7Gt39QLj9FtCYqlBmdouZdZpZZ3d3dxrVGNvbQ0+5mbY5v2kGXYdOEI0N5+TzREQmKq+T2WY2HfgR8AV3f9ef1x48Ci7lTZHc/X53b3f39oaGhkmXZSRQTMlVj+KsGQzFna5DJ3PyeSIiE5VOoNgHLEx4vyDclzKNmVUC9cCRsc41syqCIPEDd/+XhDQHzawpTNMEHEq3MpORyzkKgLam4FYeGn4SkUKXTqDYBLSa2WIzqyaYnO5IStMB3BRu3wA8HfYGOoA14aqoxUArsDGcv3gA2O7uXx8jr5uAH59ppSbi9BxFDlY9ASyeO50pVRW6QltECt64D19w95iZ3QY8CUSAB919q5ndBXS6ewdBo/+wmXUBPQTBhDDdemAbwUqnW909bmbvAz4NvGRmm8OP+kt3fwL4W2C9md0MvAF8MpMVHs3poaccXEcBEKkwls3XFdoiUvjSekpP2IA/kbTvzoTtAeATo5y7FlibtO+XgI2S/giwMp1yZVL/UJyaygoqKlIWKyvamur4yctv4e4EnSwRkcKjK7NDA9F4zoadRrQ1zeBY3xAHjg/k9HNFRM6EAkWofyies2GnEW1n1QPw8r7jOf1cEZEzoUAR6h8azkuPIlJhvKRAISIFTIEi1B+N5+waihFTqyO0zpvOlr0KFCJSuBQoQoOxOFOrcv/ruHhBPS/tO06wmlhEpPAoUIT68zCZDXDRgpn0nIqyXxPaIlKgFChC/UPxnF2Vneji5mBC+6W9x3L+2SIi6VCgCPUP5X6OAuC8pjqqIqZ5ChEpWAoUoYFofnoUNZURljbWaeWTiBQsBYpQ/1B+5iggmNDeslcT2iJSmBQoQvmaowC4qHkmx/uH2NPTn5fPFxEZiwIFMDzsDAwNU5OnQHHxgmBCe8s+TWiLSOFRoACi8eApc1PycB0FwNLGOqojFbykCW0RKUAKFMDgUBgocnyvpxHVlRWc31THi1oiKyIFSIGC4KpsgJo89SgAli+cyZa9x4nF9QxtESksChTAYCxonGvy1KMAuOzsWfRF4+w4eCJvZRARSUWBgoQeRWX+fh2XLZoFwAtvavhJRAqLAgUwMDTSo8jfr2PBrKnMnV7Db944mrcyiIikklbLaGarzGyHmXWZ2e0pjteY2SPh8Q1m1pJw7I5w/w4zuzZh/4NmdsjMXk7K68tmts/MNoevD0+8euk5PfSUp+WxAGbGZYtm8sKbChQiUljGDRRmFgHuBa4D2oAbzawtKdnNwFF3XwLcA9wdntsGrAEuAFYB94X5AXwv3JfKPe6+PHw9MUqajCmEoScI5ileP9LHkZODeS2HiEiidFrGFUCXu+9y9yiwDlidlGY18FC4/Siw0sws3L/O3QfdfTfQFeaHu/8C6MlAHSbt7cnsPAeKcJ7iN5qnEJECkk7L2AzsSXi/N9yXMo27x4DjwJw0z03lNjPbEg5PzUqVwMxuMbNOM+vs7u5OI8vRDQ7lf9UTBFdoV1aYhp9EpKAU4mT2t4FzgeXAAeBrqRK5+/3u3u7u7Q0NDZP6wEK4jgJgSlWEtrNmKFCISEFJp2XcByxMeL8g3JcyjZlVAvXAkTTPfQd3P+jucXcfBv6ecKgqmwpl6AmC4acX9+jCOxEpHOm0jJuAVjNbbGbVBJPTHUlpOoCbwu0bgKc9uGd2B7AmXBW1GGgFNo71YWbWlPD2Y8DLo6XNlEK44G7EpYtm0j8UZ/sBXXgnIoVh3EARzjncBjwJbAfWu/tWM7vLzD4aJnsAmGNmXcAXgdvDc7cC64FtwE+BW909DmBmPwR+BSwzs71mdnOY11fN7CUz2wJ8APhvGarrqAaHCmPoCWDF4tkAbNh9JM8lEREJVKaTKFyi+kTSvjsTtgeAT4xy7lpgbYr9N46S/tPplCmTCmnoqal+Kotm17Jhdw+fvfqcfBdHRKQgJ7NzbiRQVEcK49dxxeLZbHq9h+FhPfFORPKvMFrGPBuMxamprCC49CP/rjhnDsf6hnSDQBEpCAoUBNdRFMKw04grRuYpdmmeQkTyr3BaxzwajOXvMaipLJxdS/PMqWx8vSAuXBeRMqdAwdtDT4XkisWz2bi7h2CVsYhI/hRW65gng7HCGnoCuOKc2Rw+GeW17pP5LoqIlLnCah3zJJijKJyhJ4ArFs8B4Fe7NPwkIvmlQEE49FQAF9slOntOLWfVT+E/ug7nuygiUuYKq3XMk0Jb9QTBg4yubm3gua7Duu+TiORVYbWOeRJMZhfW0BPA1Uvn0jsQY8u+4/kuioiUMQUKCnMyG+C3zp2LGTz7qoafRCR/Cq91zINCu45ixKxp1VzUXM+zOyf3YCYRkclQoCC4e2wh9igArm6dy2/2HKN3YCjfRRGRMlWYrWOOFerQE8DVrQ3Eh51fvabbeYhIfhRm65hjQaAovKEnCJ54V1sd0fCTiOSNAgWFeR3FiOrKCq46Zw4/39Gt23mISF4UZuuYQ/FhZyjuBTv0BLDy/Eb2Hu3n1YO6nYeI5F7hto45Ei2g52WPZuX58wD4t+0H81wSESlHaQUKM1tlZjvMrMvMbk9xvMbMHgmPbzCzloRjd4T7d5jZtQn7HzSzQ2b2clJes83sKTPbGf6cNfHqjW8wFj4vu4B7FI0zpnDJgnqe2qZAISK5N27raGYR4F7gOqANuNHM2pKS3QwcdfclwD3A3eG5bcAa4AJgFXBfmB/A98J9yW4HfuburcDPwvdZc/p52QU6RzHig+c3snnPMQ6dGMh3UUSkzKTTOq4Autx9l7tHgXXA6qQ0q4GHwu1HgZUWPFd0NbDO3QfdfTfQFeaHu/8CSHVr1MS8HgKuP4P6nLHBocIfegL4YFsjAE9vP5TnkohIuUknUDQDexLe7w33pUzj7jHgODAnzXOTNbr7gXD7LaAxVSIzu8XMOs2ss7t74ktHi2HoCeC8+XU0z5yqeQoRybmCbh09WA+ack2ou9/v7u3u3t7Q0DDhzzg99FTggcLM+FBbI8/uPMypwVi+iyMiZSSd1nEfsDDh/YJwX8o0ZlYJ1ANH0jw32UEzawrzagKyOtZyukdRgPd6SnbdhfMZjA2rVyEiOZVOoNgEtJrZYjOrJpic7khK0wHcFG7fADwd9gY6gDXhqqjFQCuwcZzPS8zrJuDHaZRxwkZ6FNWRwu5RAFzeMpvGGTU8vuXA+IlFRDJk3NYxnHO4DXgS2A6sd/etZnaXmX00TPYAMMfMuoAvEq5UcvetwHpgG/BT4FZ3jwOY2Q+BXwHLzGyvmd0c5vW3wIfMbCfwwfB91oxcR1Fd4ENPABUVxocvauLfd3TrJoEikjOV6SRy9yeAJ5L23ZmwPQB8YpRz1wJrU+y/cZT0R4CV6ZQrE4biwRRIoc9RjPi9S87iu8+9zlNbD/L771mQ7+KISBkojtYxi0Z6FFVFMPQEcOnCmTTPnMrjW/bnuygiUiaKo3XMoqH4SKCwPJckPWbGRy5u4tmdhzl6Kprv4ohIGSj7QFFMcxQjPrr8LGLDzmPqVYhIDhRP65gl0XjxrHoaccFZ9bQ1zWB9557xE4uITFLxtI5ZUow9CoBPti/g5X29bNvfm++iiEiJK67WMQvenqMorl/F6uXNVEcq+Ofn1asQkewqrtYxC4q1RzFrWjUfamvkX3+z73QdRESyobhaxywY6VFUVhTHqqdEn2hfwNG+IT2nQkSyquwDxWB8mOrKCoK7oheXq1sbWDBrKg/96vV8F0VESljZB4qhmBfViqdEkQrjpqta2Li7h637j+e7OCJSooqzhcygaDxedPMTiT7ZvpCpVREe+o/X810UESlRxdtCZshQzIvmquxU6mur+Phlzfzr5v306EptEcmCsg8U0XCOoph95r0tRGPDfP/Xb+S7KCJSgoq7hcyAaHy46K6hSNbaWMc1583ju8/t1tPvRCTjiruFzIBobLhoJ7MT3fqBJRztG+KHG9/Md1FEpMQUfws5SUMlMPQE8J6zZ/Hec+dw/y92MTAUz3dxRKSEFH8LOUml0qMAuO2aJRw6MaibBYpIRqXVQprZKjPbYWZdZnZ7iuM1ZvZIeHyDmbUkHLsj3L/DzK4dL08z+56Z7TazzeFr+eSqOLahEpijGHHVOXNY0TKbb/6sS3MVIpIx47aQZhYB7gWuA9qAG82sLSnZzcBRd18C3APcHZ7bBqwBLgBWAfeZWSSNPP/C3ZeHr82TquE4onEviaEnCB5q9KXrzuPwyUH+4dnd+S6OiJSIdFrIFUCXu+9y9yiwDlidlGY18FC4/Siw0oJ7YqwG1rn7oLvvBrrC/NLJMyeisdLpUUAwV7Hqgvnc/4vXOHxyMN/FEZESkE4L2QwkDnrvDfelTOPuMeA4MGeMc8fLc62ZbTGze8ysJo0yTthQfJiaEulRjPiLVcsYiA3z9adezXdRRKQEFGILeQdwHnA5MBv4UqpEZnaLmXWaWWd3d/eEPyzoURTvldmpnNswnT+46mx+uPFNNu85lu/iiEiRSydQ7AMWJrxfEO5LmcbMKoF64MgY546ap7sf8MAg8F2CYap3cff73b3d3dsbGhrSqEZqpbI8NtkXP7SUeXU1/Pf/+xKxuJ5XISITl04LuQloNbPFZlZNMDndkZSmA7gp3L4BeNrdPdy/JlwVtRhoBTaOlaeZNYU/DbgeeHkyFRxPqc1RjKibUsX/+EgbW/f38rBu7SEik1A5XgJ3j5nZbcCTQAR40N23mtldQKe7dwAPAA+bWRfQQ9DwE6ZbD2wDYsCt7h4HSJVn+JE/MLMGwIDNwOczV913K4V7PY3mdy9q4p+X7uWrP93B+5fNY/HcafkukogUIQv+8C9u7e3t3tnZOaFzl/3VT/jMe1u448PnZ7hUheHA8X5W/Z9naZk7jUc/f1VJ9p5EZGLM7Hl3bx8vXdm3GqU6RzGiqX4qaz92IS/uOcbfPd2V7+KISBEq3RYyDbH4MMNOyf+V/ZGLz+Ljlzbzrad38uzOia8QE5HyVNot5DiG4sGwWyn3KEb89fUXsmTedP7kh79hT09fvosjIkWk9FvIMURjwbLRUu9RAEyrqeT+T7czPOzc8vDzuheUiKSt9FvIMUTD6wuqS+yCu9G0zJ3GN2+8lFcPnuDz33/+dKAUERmLAgXlMfQ04v3L5vE3H7+IZ3ce5ovrNzM8XPyr3kQku8a9jqKUDZXR0FOiT7Yv5OipKH/zk1eorqzgq79/MZVl9jsQkfSVdaAoxx7FiM/9zrlEY8N87alX6RuM840bl1NTGcl3sUSkAJVfC5mgnCazU/mTla3c+ZE2frr1Lf7wu5s4eiqa7yKJSAEqzxYyVM49ihF/9L7FfP2Tl9D5+lGuv+85Xj14It9FEpECU74tJG/PUZTKM7Mn6uOXLWDd566kLxrn+nuf45FNb1IKt3YRkcwo6xZSPYq3XbZoFo/d9j6WL5zIUIvSAAAIpklEQVTJl370Ep97+Hm6T+gJeSJS5oFiKF7ecxTJ5tdP4fs3X8Ff/e75/HxHN9d87ed877ndep6FSJkr6xYyqqGnd6moMD579Tn85AtXs3zhTL782DZ+95u/5KcvH9A1FyJlqqxbyOjpez2Vx5XZZ+Lchun84x+t4Dv/+TKi8WE+//0X+Mi3fsnjW/af7omJSHko7+soTvcodP1AKmbGqgub+OD5jXS8uJ9vPd3Fbf/0G+bV1XDjikV86vKFnDVzar6LKSJZVtaB4vQchXoUY6qMVPDxyxawenkzP99xiId//Qbf+NlOvvGznbzn7Fl85OImVl04n6Z6BQ2RUlTWgUJzFGcmUmGsPL+Rlec38uaRPh7bsp/HXtzPVx7bxlce20brvOm8r3UuV7fO5bJFs5hZW53vIotIBqQVKMxsFfANgudb/4O7/23S8RrgH4H3AEeAT7n76+GxO4CbgTjwp+7+5Fh5mtliYB0wB3ge+LS7Z+WS4bd7FAoUZ2rRnFpu/cASbv3AEroOneCZV7r5xc5u/mnDm3z3udcBaJlTyyULZ3LxgpksbZzOuQ3Taaqfgpl6cCLFZNxAYWYR4F7gQ8BeYJOZdbj7toRkNwNH3X2Jma0B7gY+ZWZtwBrgAuAs4N/MbGl4zmh53g3c4+7rzOw7Yd7fzkRlkw2GPYoaBYpJWTKvjiXz6vjj3z6HgaE4L7x5lM17jvHinmNs2NXDjzfvP522tjrCuQ3TWTS7lqb6KTTNnErzzCk01U9l3owaZtVWM6VKc0YihSSdHsUKoMvddwGY2TpgNZAYKFYDXw63HwX+zoI/G1cD69x9ENhtZl1hfqTK08y2A9cA/ylM81CYb1YCRV80RmWFaegpg6ZURXjvuXN577lzT+87dGKA1w6doqv7JK8dOslr3SfZfqCXf9t+8HSwTjS1KsKs2ipmTatm9rRqZkytYlp1hNrqSqbVhD+rI9TWVDKtupIpVRVURYJXdaVRHYlQVWnB+0gF1ZUjx41IhVFhhhlE7O1t9XJERpdOoGgG9iS83wtcMVoad4+Z2XGCoaNm4NdJ5zaH26nynAMcc/dYivQZd2owztTqiBqJLJtXN4V5dVO46tw579jv7hztG2L/sX4OHB/g8MlBek5FOdYXpefUUPCzL8q+Y/30R+OcGoxxKhonnoXrOSoMKsyoqLDT25EwiFRUjGwHx0aM/LMxLMW+kffv/rdlY+TxjmNJebwjp6TPKWfl/v/3f3/sIlYsnp3VzyjayWwzuwW4BWDRokUTyuP8pjo+fGFTJoslZ8DMmB32Gi5srk/rHHcnGh+mbzDOqWiMU4NxBmNxhuLDDMaGGYo7Q7FhhuLDROPDREf2hdvD7gw7wc/hYDvujrsTD98nbgfpg1cwpeVhOXjHz+BI0rF3lPudaRIT+Ok0nnxolLz8XfvKln4JTKvJ/lBtOoFiH7Aw4f2CcF+qNHvNrBKoJ5jUHuvcVPuPADPNrDLsVaT6LADc/X7gfoD29vYJ/XP51OWL+NTlEwsykh9mRk1lhJrKCLOmaVWVSC6kMzi/CWg1s8VmVk0wOd2RlKYDuCncvgF42oM/ezqANWZWE65magU2jpZneM4zYR6Eef544tUTEZHJGrdHEc453AY8SbCU9UF332pmdwGd7t4BPAA8HE5W9xA0/ITp1hNMfMeAW909DpAqz/AjvwSsM7P/BfwmzFtERPLESuG5A+3t7d7Z2ZnvYoiIFBUze97d28dLp3WhIiIyJgUKEREZkwKFiIiMSYFCRETGpEAhIiJjKolVT2bWDbwxwdPnAoczWJxioDqXB9W5PEymzme7e8N4iUoiUEyGmXWmszyslKjO5UF1Lg+5qLOGnkREZEwKFCIiMiYFivDGgmVGdS4PqnN5yHqdy36OQkRExqYehYiIjKlsAoWZrTKzHWbWZWa3pzheY2aPhMc3mFlL7kuZWWnU+TNm1m1mm8PXZ/NRzkwxswfN7JCZvTzKcTOzb4a/jy1mdlmuy5hpadT5/WZ2POE7vjPXZcw0M1toZs+Y2TYz22pmf5YiTcl812nWN7vfs4dP9yrlF8GtzF8DzgGqgReBtqQ0/xX4Tri9Bngk3+XOQZ0/A/xdvsuawTr/NnAZ8PIoxz8M/ITgCaJXAhvyXeYc1Pn9wOP5LmeG69wEXBZu1wGvpvi3XTLfdZr1zer3XC49ihVAl7vvcvcosA5YnZRmNfBQuP0osNKK+2G86dS5pLj7LwiehzKa1cA/euDXBE9TLOpn4aZR55Lj7gfc/YVw+wSwHWhOSlYy33Wa9c2qcgkUzcCehPd7efcv+nQaDx7DehyYk5PSZUc6dQb4/bBr/qiZLUxxvJSk+zspNVeZ2Ytm9hMzuyDfhcmkcIj4UmBD0qGS/K7HqC9k8Xsul0AhqT0GtLj7xcBTvN2jktLxAsFtGi4BvgX8a57LkzFmNh34EfAFd+/Nd3mybZz6ZvV7LpdAsQ9I/Gt5QbgvZRozqwTqgSM5KV12jFtndz/i7oPh238A3pOjsuVLOv8OSoq797r7yXD7CaDKzObmuViTZmZVBI3mD9z9X1IkKanverz6Zvt7LpdAsQloNbPFZlZNMFndkZSmA7gp3L4BeNrDWaIiNW6dk8ZsP0ow9lnKOoA/CFfEXAkcd/cD+S5UNpnZ/JG5NjNbQfB/vpj/ACKszwPAdnf/+ijJSua7Tqe+2f6eKzOVUSFz95iZ3QY8SbAa6EF332pmdwGd7t5B8EU8bGZdBJODa/JX4slLs85/amYfBWIEdf5M3gqcAWb2Q4LVH3PNbC/wP4EqAHf/DvAEwWqYLqAP+MP8lDRz0qjzDcB/MbMY0A+sKfI/gAB+C/g08JKZbQ73/SWwCEryu06nvln9nnVltoiIjKlchp5ERGSCFChERGRMChQiIjImBQoRERmTAoWIiIxJgUJERMakQCEiImNSoBARkTH9f1gspyx5meS8AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Define a discrete-time power spectrum with sample time Ts = 1e-3, bandwidth 10 rad/sec, and \n",
    "# standard deviation 3 (when filtered by white noise)\n",
    "\n",
    "wu = 10 # bandwidth of the force disturbance\n",
    "std_du = 3\n",
    "Ts = 1e-3\n",
    "Hu = control.TransferFunction([1], [1/wu, 1])\n",
    "Hu = Hu * Hu\n",
    "Hud = control.matlab.c2d(Hu, Ts )\n",
    "t_imp = np.arange(2500)*Ts\n",
    "t, y = control.impulse_response(Hud,t_imp)\n",
    "y = y[0]\n",
    "plt.plot(t,y)\n",
    "std_tmp = np.sqrt(np.sum(y**2))#np.sqrt(trapz(y**2,t))\n",
    "Hu = Hu/(std_tmp)*std_du\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2.6940445662402728"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJztvXm4LGdV7/99u2voufd85inJyTySk4EEAiSRMRqvouIl/EDEIAiC+oMfwiPqVRQBRbwqGhmFXBAjCBchQIAQIhByMpHpJOfkzOOed89dQ7+/P6re6urqqu7q3VU9vp/nyZN99u7dXbuHVau+a63vIpRScDgcDmd0iPT7ADgcDocTLDywczgczojBAzuHw+GMGDywczgczojBAzuHw+GMGDywczgczojBAzuHw+GMGDywczgczojBAzuHw+GMGEI/HnRmZobu3LmzHw/N4XA4Q8tDDz20SCmdbXe7vgT2nTt3Yu/evf14aA6HwxlaCCFH/NyOSzEcDoczYvDAzuFwOCMGD+wcDoczYvDAzuFwOCMGD+wcDoczYgQS2AkhE4SQuwgh+wghTxNCnh/E/XI4HA6nc4Jqd/wYgLsppa8mhEgAEgHdL4fD4XA6pOuMnRCSBXADgE8CAKVUoZSudnu/HA4nPE6ulvHtJ0/3+zA4IRGEFLMLwAKATxNCHiGEfIIQknTeiBByOyFkLyFk78LCQgAPy+Fw1subPrsXt3/uIeQqar8PhRMCQQR2AcDzAHycUnoFgCKA9zhvRCm9g1K6h1K6Z3a27UQsh8MJkaKiAQC+9ujJPh8JJwyCCOzHARynlD5g/vsuGIGew+EMKLvnUgCA7+2b7/OR9J6Kqvf7EEKn68BOKT0N4Bgh5DzzWzcBeKrb++VwOOFR1WoAgKdO5vp8JL3l+EoJF77/bnz54eP9PpRQCaqP/e0A7iSE/AzA5QD+IqD75XA4IZArG9r66VwFVW30M1jGvlN51Cjw6f8+3O9DCZVA2h0ppY8C2BPEfXE4nPDJVTTr6/lcFdumxqNDeT5f7fch9AQ+ecrhjCFrZRXbpuIAgPl8pc9H0zsWzMBOSJ8PJGR4YOdwxgxKKXJlFefOpQEAp9fGI4ullOKj9zwLAKiqtT4fTbjwwM7hjBllVYdWozh3oxHYjy6X+nxEveHpU3nr6/yI9+/zwM7hjBm5sqGvb5tMYCop4djKuAT2egdQvqq1uOXwwwM7hzNmFMyglpSjmEpKWC4ofT6i3vDosVXIQgRvffHZKFQ1UEr7fUihwQM7hzNmlBWjvTEhCZhJSVgojIfG/tSpHC7fNoFMXASlhiQ1qvDAzuGMGSygJaQoNmZiOL02Hl0xp1bL2DIRR1I2urwLIyzH8MDO4YwZJdMnJiZGsTEbx3y+glptdGUJANBrFGfyVWyaiCElRwEAxSrP2DkczohQl2Ki2JiRoeoUy6XR1tnXyir0GsVMSkZSMjL2Is/YORzOqFCyB/asMaQ06nKMXX5iUgwP7BwOZ2RgQS4uRrExGwMA3PXQceg1imMj2tPOrlJioi2wKzywczicEYEFubgUxZYJI2P/zI8O49P/fQgv/ND3R9LxsWI7mTGNvcA1dg6HMyowKSYuRjGblpGNiwCA+/YvAgC++cSpvh1bWFhXKVyK4XA4o0hZ1SFFIxCixsf/ddfuQDRCrM6Yw0ujJ8eUbSezBC+ecjicUaOsaIhLUevfmyZi0GsUB+YLAIAzI1hIZRl7TIwiKfF2Rw6HM2KUFB0JW2DfzDpjckZAd1tw/YkfHsSH7t7XmwMMgYqtK0aIRhATI7x4yuFwRoeSqiMu1gP7dEpq+LlzInMhX8Wf/9fT+Md7n+vJ8YVByVYwBoCULPDJUw6HMzpUFL1BiplMNAZ2p/Z8/4EF6+thNc6ya+wAkJQFlHhg5ww6X/zpUfyfB472+zA4Q4BTiplK1gN7NEKatOcfmt0ywPAaZ9k1dsAwQOPtjpyB5thyCe/58uN471ce7/ehcIaAkqpbAQ5AQ5DflI1B0WsNC65/8tyS9fXSkFr8VlQdhACyYIS8TExwrSWMCoEFdkJIlBDyCCHk60HdJ8cfr/vkA9bXw3qpzOkdZUVrCObEtgB0szmwxLL2Wo3idK6C881tS0vF4QzsZcWoK7C/dSopYWVI/xY/BJmxvwPA0wHeH8cn9r7jlZJ3FlKoavib7zyLtRa34Yw+hYqGlCy6/myTaTHAdPZ8RUONArs3GIF9MT+c3u1lR8F4MilhZYSNzwIJ7ISQrQBeBeATQdwfpzNiYv1lPLla9rzdVx4+jr/77n78w70HenFYnAElX9WQjgmuP9tktj7mK0ZgZ8Hv3LkUAGBxSJdylJVG+WkqIWG5qIzsFW5QGfvfAng3gNFe/T2giJEIrtg+AQA4k/MeLvnbe/a3vQ1ntKnVKAotAvvmCTNjVxoD+zlDHthLio6kXA/sqZiA2ghvUeo6sBNCbgEwTyl9qM3tbieE7CWE7F1YWGh1U04H1GoUBUXD9qkEAGDVQ2Y5k6tY+ujyCGuLnNYUFQ2Uoimwv+XFZwMAds8ZkkvBzNjZ+2lDNoa0LGBxSIunJVVHXKr/zaO+RSmIjP16AL9ACDkM4IsAbiSEfN55I0rpHZTSPZTSPbOzswE8LAcACuYHddukEdi9dMPnFgr13xnRNzOnPey1T8caNfZ3v+w8PPr+n8OMOay0UKjit/51L/YeWQZgSBczaXlo96OWqpplJQDA+ro0oi2PXQd2SukfUkq3Ukp3AngNgO9RSm/r+sg4vmBa6KaJGCLEO2NfMIteu+dSI21+xGkNe7+k5MaMnRCCiYSElJnJ/+jAIr7z1Bn8w/eNadOplIS5tIzv75vHW+98aOjeQ0bvPs/YOUMCu2TOxkVk4yJWy+4Ze65sBPwtk3Hrdzj9h1KKW//hv/Gdp86E/jg73/Nf+KP/fAJAsxTDYAH/ZyfWGr6fiYk4ey6FkqLjG4+fxuOOnw86ZbVxKCvNA7t/KKX3UkpvCfI+Oa3Jm0MW6ZiIyYTk2e6YY5l9No78iL6Zh5FHjq3isWOr+K1/3Rvq47CW2AcOGdKKV2Bn+0APLhSbfsbqOMBwFeBrNYqFfBWTibr8xK5MRjXJ4Rn7kMMurdMxARMJEaseGnuuokISIphJSShWtZFt8xo27nvWaCS4cFMm1Mc5tFho+LdTY2dEIqTpe29+0VkA6sNLwHDtSD28VEShquGizVnre+zvz1dHc6bD/bTNGRrYWHQmJmAiIXlmUrmyhkxMREqut3nZNUdO76lqutWCWgv5RPvkicZ1d14ZOwC84JwZ3H9gEW+/8RyUFB2/d/O5AIANadm6zekhythZJ4/9xMQkp1HN2Pkne8ipZ+wiJhIinjmdd71drqIiExPqRaOK1hTYi1UN0QhpGOQYRyil+P0vPYbbrt2OK3dMhfY4P7Z5sITdT73vTOP7wlk8tfP5N11jjOBLje+Di7dk8bKLNuBbT54ZKimG1ZfsJzP2dW5EAzuXYoYcuxRjaOzuUky+oiEdF60CEvOntnP9X30PP/fRH4R3sEPCfx9YwlceOYFf/viPQ30cNk9w9mwy9G0+Ts082eZqzRnUAaOT5J9ftwfXnT09VFIMk1sy8br8JAsRiFFifX5GDR7Yh5x8RUU0QhAXo5iIiygpeoMzHyNXNjL2VoF9taTi2HJ57PX3R46uWF+XQtqys1ZW8cBBo5B5/TkzoT0O4+hSEVfumLT+7aal+2VjJoYzueHpZ8+V68kPgxCClCxYzQejBg/sQ84TJ3PQa9ToQzZ9td1MvnIVFZm4aE3ftbr0b2UkNg7YW/3mQwpg7/r3x/Bve49hLi1jMiGhpOjWMumgUfUaioqOCzalA7m/DdkYzuQqoR1v0LhJMca/Rd7uyBlMnj5VL4qxdi63wMyKp8zhruzI2DW9bvOzNKTThUFx37MLlgY9H5Kb4bfNvvW4FLU8TMLS2ZncsHM6Gcj9bczEoNXo0Fj45qsaYmIEstAoL6VjApdiOINJsarhVZduAlBfceams7PiaV2KaXxD29/g45yxH1wooKrVcJ7pP74Qsk3tLZdusq6iwlquzOSGiYSEndMJ3HzBhq7uj21c8mqtHTRyZdW1vdMI7KP5XueBfYhZK6soKTou2WL0506YGbvzA1dWdChaDZm4aHW8OLND+zaZUfapbsfvf+kxAMBvvXAXgOb9n0FAKYUYJXjphRvwuzfttnxLwiqg2gvs33zHDfj4bc/r6v7YFUbRpU4ziOQrGjIu7Z3pmDiyGTtvdxxiDswbQye7TUtVlkktFxuzkPm80cEwl5atjN0pxbACEzA8mVgYHFkyukeuPWsaQDhZdK6sQdUprt41BVmIWvJYJSQpJlepa8xu3S6dwtpkh2UZdK7SKmMfjr+hU3hgH2IOzBu9ycxq1U2KqdUo3vw5w1F5Qybm2RVjz9i9jMTGgbgYxY3P21APXiFkpQdMp81Zc+AnFnJgZ8Er4zFt2ims/jAsGXuurCKbkJq+n5KF0LuR+gWXYoaY/WcKiIkRbJk0JupiYhQxMdKQcT9+Yg37zKGlrZNxbymmXA/mo5rF+IFtF5LMPucwpJinThpdNyywswXLFTWcPTV2KSYIvOo0g0q+4r5YJCEJQ3Ny6hQe2IeY07kKNmZiiNp6kiOE4KEj9T7sQ4v1wZQtk3HIQgQR4iLFVOyBfTwzdud2oYQkhJKxq7rRJsj8YWTzZOs2fxAEdqO4IAjzaiYMjMaB5r89KUWhaDWo+ugtfuOBfYhZc7nELCk6Hj66ilqNYu/hZbzz3x4FAPzoPTdCFowt7XEx2izFmBr7RGJ0C0rtKKl6w3ahpBQNJWNnvdPM3oHtrA07Y29lI9AJCVY8HSCNnVKKO+57DseWS00/88zY5eE6QXUCD+xDzFpZxUS8MRP5jet3AgA+86PD+PC3nrG+bzdAiksCymrjhzJXUUGI0aM8qv4Z7WCZbUo2ntO41HwCDOpxYmIEYtT4+LH+6rAy9lxZRVyMQhKC+bgnRO/p5X7xzJk8/uIb+/DWOx9u+L5eo6hqNatAbSc5ZJJSJ4x8YB/FyyzGWlm1WhwZt99gWKz+yw8P4pFjq66/l5KjTa11ubKKtCwgExcbZJlecHylhHff9VhoxUO/OLVoQ4oJJ2O3yyIsY6+GlLEbU8fB9UkI0QhkIRJa3/16eOqkMahnXwEJ1GtJCZduIJaxh+3T0w9GOrA/cWINu9/3Tdz7zHy/DyUUVksqso6MfVPWyMxPrVWgaEag+MtfuqThNimXwYy1smE5kOmgBays6LjzgSNdnzz/6u5n8KW9x/G9ff19nZhEwiSL8DJ2zdrgA9i6YkLL2LWm90m3JGVhoPaFstfO+XqxWpJrYPeYwh4FRjawP3Uyh7/7ruF1/cP9i30+muCp1ShylWYpxskbr9+FX796e8P3UrLQ5JFxcrWCzdk4MjHRd/H0iw8exfu+8gS+/rOTnR28B2FJEX6pmB9w1uudkKKhjPkXqpq1wQeod8WEmrEHVDhlxMXoQGXs9uBsN7Fj33ezorZqBQP0dwTFyAb2N332QcuPo98BIwxyFRWUNlqRumF39GO4TdwdXS5h21Sio6GNJXOBwf4zhTa3bI1eMwIau8LoF87L9kSYGXvMJWMPcUCp3fukU5JydKAydvvrZG9hLJm1JLelMkmru2f0AvvIDiidtPlFh+XQ10+YAdNMSm762affcBVOrJbxcxduwIZMrOnnabkxeOs1ivl8BZsnYqDUyCgpNRwjW6GZ7n7r3abzwMElbJtKQDPb//rdjcOCAyu0xUUhlMv0QkXDTKq+P1SIEERIuFLMObPBftSNHvDBCYj2k+JqSbHktFILKcayRhigE1RQdJ2xE0K2EUK+Twh5ihDyJCHkHUEcWLeI0XpQGhYXuk5g2fJ0qnmi7iXnz+G2a3e4BnXAKA7apZilQhU1algOpGMC9Br1lamWzQ/2lx8+0fHx12oUv3bHT/C7X3gEiqnR97sbh2XsLIM2MvZwiqes8wYwvMFjYjTk4mkIGfsAadPlhsBelxKd8pqdxAhn7EFIMRqAP6CUXgjgWgC/Qwi5MID7XTcVVYeqU7zrZefh5y/bPHI2tJpew6/+s7HdZzrZnLG3I2UGdqZFMmva2XTM6tbw0xlj77t+/PhaRws62GPuPbKCNXPq1T792g8qPZJiDO+Sxgw6JkZDydgppeaSlWADe0ISBqqP3f462QO78yrMDpNieMbuAqX0FKX0YfPrPICnAWzp9n67gb2wkwkJ00nJym5HhSO2IYy5zDoCuyxCr1Ery2EBdTIhWm1xdlMwL+yB6Of//n7864+P+D6Go7a/4eRqGUD/pZiyI7uLS1FUtRr0ABdKUNo43cqQhUjbAaXnFgr4628/09HxFBUdNYpA2x0Bowd8UDN2e1JSatHuGOd97P4ghOwEcAWAB4K8305huyQnEyJmUhLyVW2kCqjMMuA1V21z1djbwYIK29Bun4SciJtbmHxkz0792b6cuR0nVuuBna1Z63X/vBNLihHqGTsQ7Ae/pBjTrc4p0JhonERa8bf37Mf//t4BPG7b8NQOdtIOvCsmJLuF9VJRdMvd1H7l10qKkYQIpGhkJP1iAgvshJAUgP8A8E5Kac7l57cTQvYSQvYuLCwE9bCuMHfDyaRkBb5RKqA+dNgI7H90y/oULxbY847e34QUtfqd/Vj3VrRaQ4DqpJ/95GpzwbXfUkxZ0Q0vHdN7x1ojGOAHvz4E1RhojYy99ePsM7dlsWEcP7CTZeAae0j1h/VSUnSrptSQsSveXTGA0fI4LPbDnRBIYCeEiDCC+p2U0i+73YZSegeldA+ldM/s7GwQD+sJy9inkhK2ThrdB8dXyqE+Zq84tlzCv+09hgs3ZSyvkU5hwZgFGdbdkJIFK7D7ydgrqo6zZuvr1jpZI3dqrYxsXIRgMzDrd/G0ouoNmV0Yo/OFqmlb4JRifGTszJb5vV953HdrJJPUAtfYZSHUPa2dUlZ1zKQkREjje5dJMW4aO2Do7Dxjd4EYPXGfBPA0pfRvuj+k7mHZ5mRCwrYpYxLz2EqzOdAw8pkfHQbgL/B6wbJFJsWwfuTEOgL7pM2EjC2p8MOB+QJ2zSTxxduvxfkb03jh7pn+Z+yq3hAArKUkAfaXWxm7U4rxkbGrtXrgP7HqL1GxpJgQNHYgvD2tnVJWdCSkqGGJYasPlRUdhNRtG5zEpSifPPXgegCvA3AjIeRR879XBnC/64ZtEJpIiNg8EUeEAMddXN+GESZ3fPhXLl33fdQzduN5Yhp7QowiHRNAiD9ZpOIIhLmK5mtH6EpRwaPHVnHp1iz27JzC3e+8AefMpfoe2CtqrWFCMe6xlKQb2HPdVDwVo6i2CZL24vJJv4G9Eo7GbvmsDIgcw07KWYfXUVkxvu81k5GU6hO0H7tnP/7+e/t7crxhE0RXzP2UUkIpvZRSern53zeCOLj1slJSkIkJEKOGg96mbBzHRkSKWSooOGs2ievOnln3fTg19rKqIyYa2nIkQpCJiVj1FdhriIkR/NfvvgBvffHZAICrPnBP29977PgqKmoNr7h4k/W9TExEvqoF2oHSKcbzYM/YjecpyIlQy0LX2e4oRNpKMYWKZk0SL/uczahn7MFr7AD6Mn360JFlvPJjP2yoA5VNGS0TExsShJKqu3bEMBJS3fPmo/c8i498+1loI2AcOJKWAstFxaqQA8bmIKfr27CykK9iLt15J4wdZ1eMM/POxkV/XTFmILxocxZvftHZ1vcv+ZNv4cmT7p0ba2UVn7z/EABg+3R9+pIFnn4u+aiYJzhGPAyN3cMbPSZG255AClUNW0z7Zb9SXC7g7UkMdtLrR8b+wW/uw1OncnhuoS79GZm5gExcaHhuyoru6hPDSMpGxm6fwbDf77AykoF9paRg0hbYJSGCnx1fw6f/+1AfjyoYgjB0SjqKp843v9/AXrFluNm4iP9160XW/X7E5gVv58+//pRlyrbBdoLqRNsPC+cJLow+53y1VVeMd6aom9udWM3I717aXFlFQopa3u9Bwcbx+9HyaA20mUkApdTM2CNGxm6TrEqK1j5jV3TrdQGAZ8/kQzry3jGSgX25qDQU9abNIP/VR4NxIewnuXL34+FiNIKYGLE6NCqORQR+DZ6qDk161tZTf8DjCsl+5STYgg0L7H4Go8LCS2MPVophyzzc+ti9H4f93nRSRlKK+g7shn1B8JZQVsbeh1ZBJlkxyUXVKfQarWvs9oxdrVltq24kZWNLlr0d+nM/8T9oN6iMZGBfKSoNCyj++OeNTNIuzwwr+YoWSCEsHROtQl5F1a29mwBrAWv9gdVrFIpea5Au7FOwx5bLrpmul4KeMaWCXmbszsdydsUwO90gXSdXSyoyMaFhTy1gdG20ytittsW4iImE5Pt5Ciuw9zNjZ4+Zs11xAsbcgXNRTFnRrLZVN1jGbs/Sf3poOYzD7ikjF9hVvYYz+aqlRQLGoNILd89Yg0vDil6jyLuMo6+HtCxYHwyntsx6lFtRcekPnk01mo45h8I0vWYN1/z5L17c8LNsordSzBd+ehSX/em3G5awGCe4+vNg+aQHGNhXSgomEs0JhiwYXjFefjvsecnGRWTioq8BMsDIqNc779CKZJ8y9qVC1eq8YlcxZdt7MRMTUFFr1tVPSdFdp04ZrCvmiRNrECLEmsvo9zavbhm5wL5aUqHXaFOBcSIhYWXIXR5Z4S2IDge777pTW/azxNnphAgAG7IyZtMyfm3PNgDNA0sLhSqqWg0f+B8X47ZrdzT8zJJielA8PbFaxh9++XEAwA+erU9BO58HKYTAvlpSMZlofv1iYgSUGrKCG/W2RaHJnbMVxapuZddBEkaPvx/s2TS7iqlPl0abJL1ym8CekAVQalhPZ+Mi3vQCY7WkX6lrUBm5wF7wKE5NJUTfLWKDCvtwB5GxZ+L1TUlObTnhwweEafD2bFAWonjwfTfjN16wEwAwn2+0DThteuRvdLETZvJSLzJ2+yCVve/e+TxI0TACu3vG3m49HtONswmxIwOu8KQYwbr/XrJoc2rNOTL2mBi1kh72s5Kit5RiWNvmqVUjsDMJd7U83LFi9AK7RzvZZFJCrqINdY9qkMMm9n7fkqI1FU+dLWBOmAafdMmG5tJG4D65Wm6YkDxjLuRw84lPSFEQ0ptLe7vvOTvZsM4KuyRFCIEkRALV2BcLilXMt9NuPR474WViIhKy/yUXRSUcKYYdbztHyqBZKCggBNg2FW+44gRg9bED9enzUpu/nxVWT62VkU3UA/tKkWfsA4XVdeDIalnh1M/gzaBS9/0IQGOP1TX2lZLaUGxOSMblaasPbdHmCOlkMiFCjBL8xTf24foPfs8qblkZe7Y5sBNCEBd7M97NAsFZs0nrKk7RDXtep1mULEQCcwallBpzCC4nNla89tJ27WZeScn/WrpiVfc0wOoGQgikaOuT3sGFQuBeMouFKiYTEqYSkpWYMD/1hBS13lu//PEfW9bU7TR2wNi4lo2LNndTnrEPFKwftSljNy9/h1lnZ5e9zpPWemBSjF6jWCk1ZpFJH0t+mXGSm35LCMFmW/H6wLzR4ngqV4EYJZhykSIAc7FFDzRbJq1smYhb27WszgrHZbscYMa+WlKh6DXXAbN6odZLitEQjRAkpWhHa+mKVQ2pEDR2oPVJ77Fjq7jxr3+AP/7ak4E+5nyuYm76qkuJRdtn/rwNaeu2ew8vQ9Vp664YM04oWq1RiuEa+2BhFRgdcgUL7MOss7ezIO0E1j2wkK+CUmDa1oNurQxrkRUyq1OvY9k1U3d9PLhoZG6f/OEhxMSoZYvrJNbjjH1zNo61sgpVrzUtsmbIQnvXRb+wYrLbcpT6QmtvKSYTE0AIsdbStdtYxTLWMKQYAJBFbxuEQ4tGHePfHjwW6GOezlWwIRNDJl6/4izYAnskQnDP778IAPCUaXOcaPH326XEbFysJ4A8sA8WXlltvetiMEyL1kOpRZbcKay4zD6A9h5/9mZvlRUWPK6MGPYT60pRwY+eW4JmDpF4keiR054V2M2ripWSUl+h1hTYg8vYWeHPbTkKC+yeGbttb2lCMvbStjvh2O2Yw6CVFMMWnMserorrZbmgYCYlIy2LTSZ27O/cZMoxR5YM479Wk6dxR2CPiRFIQoQXTwcNr8m++sq34T0TFy0XxiCkGOM+DpsdInYphmU4rUbpW21/B4Artk9YXy+XVOw7bWRPn3j9Hs/7jEtCT6SYihmMNk0YAWC5qFgnFOcViBSgxp5v4dvSrhi5Vlat5MQy4GpzEmxVBwmCVh7yrFDeyfIVP6yY7aKZuGDVnJx/Z1IWkJCi1nu7naUAIxsXQQjBZELE8pCv0xy9wF7VIAsRqweZwTLIfppMdYtzJ2c3ZJwZe8olY28hxVhdMR5B47XX7MDnf/MaTCRELBereG6hgKmkhEu3TrjeHgDiYsRaZRYmrPNkc9bI2JcLiueJSvbhuugXqxVXdutjb5Ox2xZSW62Gba4+7UXFMDCuZtyPlxXKK2otsKuwiqqjrOqYTEpIx0SUVR2qXkO+qhlr7myf+bm0bMvY/Ukx7IooLkbx7w8d76tvUbeMXGAvVNwnM5k040eK+ecfPIdP3T94hmElVYcYJU0nrfXAsr+DCyxjd9HYWxVPq0YxT/Y4FkmI4AW7Z0xDMQ0LeaWtK2VCElBSw5fKKprxPG6dNAL78ZW6/YHTCTDIdseCR8eW8bitM/a87X3NumrO5JvXC9ophS3FtDjp2T2BlgOa+GYFzYmEaHWG5SuaabPR3AXHkpZW0qVTimG/CwD/eO+BQI67H4xcYM9X3AcyxGgECSnqS4r5y2/uw//6+lNhHF5XlKpaS426E9ibmF2u2qch/fiAGG103gsMGBmze2Gt3Ojf40a8g8GbbqioOmQhii2TxhKWYyslmxQTXvG0vjS8+TWUhdbtjiWl3rbIZLN2jQD1jD0kKUaIuPbd12oUh5dKOMssoDs70Y4sFdc12LRi24yWtl2BF1w+86fW6ie9VgmFU4oBgM++8WoAwKpHL/taScXew4PtJzNygb1Q1ZqmThlGkPH/hhqkZb2A8eE0AXitAAAgAElEQVQOSi/NWBl7ARMJscFpse613aIrRtEsv5BWMOuC1ZJq9Qh7ERejvZFiNMO8TIxGsDETw4mVsqcUE2jGXjWuFFgQt8Mydq+TiN1+lmWU7Vp3S4r3iSQIZCEKxUVDXyxUoWg1XLI1CwDYd7pusFVSNLzow/fi4j/+VsdLVVhgn0iI9QnTsoZ8RW36zLPFL4D7QBzDbsbGEo90TMT5G9NWK6yTt33hYbz6n37cMAU7aIxeYPfI2AE2lNM6Y7cXe9ppmL2m1GbYohNYdlKjza6XVsbeIqsqVnUkfAQMNuG6WlbbZuy96mNXtJoVXLdMxnFitVxfeuyqsQdzTK2mIH1l7ObzbbXutpE4ih4F4aDwKiwfN6eNr9k1DQD4f//9Mfzs+CqAxkx6qdhZYGRSjJGxM2lVbZCpGK+9pu5F5JXoOcnaPJgytj55J0+bbZQPHBzcrD2cV7yP5Coqtk0lXH/mtPR0wx7M81UNc4EeXXeUqq2XBnRCTIxa2ehMsvFSNSaY4/2tpBjFnwcJy9jXSqrl4OhFr6SYqlazagNbJuJ48PCKdRJzXoW4tTueWivjyRM53Hzhho4et1jVPa9yWmXsrLWRdUPFpShiYqTtEE2phfQTBF5SzAlzDaW9M+rpUzn8+LmlBsuIQkXDXLrp1z1hvj4zKRmshT9vBvadM42f+UiE4AfvenFHV1v2wJ6OCVbLppOZlIzFgmLJmIPI6GXsLWxtMzZHQy/sP+/HEoFW2HXWIGBdFs6MPRIhSIjRlhl7qdp6lyQjHRNxJl+Botd8STGKVgt972lV1a0C9EWbszixWsa+03mIUeIqxTiD7W9+Zi/e9K978czpzjbtGBa67s9Zq4zdTVKZSkhYatOS14uM3SnFVFQdb//CIwCMlZQ3X2CkRk+cyOEvv7kPf/e9ekGyVdeVGwv5KqIRgumkLWP3kGIAYMd0Ers3+D9z2OtX6Raxgr1GhxdHPLATQl5OCHmGEHKAEPKeIO7Ti2PLJVz1gXvw4+eWXH9eqGpIe0oxYtviqT2j7/SNFzblNot5OyVr9rLbWx0ZhtGU999fqPrT2DNxwcqu/EgxQPhWsPaM/fpzjKXgX3nkBCYSUlMx2K14yiYaO12hVlQ0zyAbjRCIUeLaFePW5jqVkrDcRsqoTweHk7HHhOaBsoePrFhfp2MiPvH6q7B7LoW7Hjre9Pud7ktdLRtL6iMR0uDi6NUw4ZcP/fKleMN1Oxtee/siGiesaH18pez680Gg68BOCIkC+AcArwBwIYBfJ4Rc2O39evG5nxzBQr6K/3zkRNPPKKXGi+yVscc7y9jLPWi964RigFIMUNceZ1zcBg1r2FYDSv5cA+2Z1EQbH3nmtBf29GlV060M+fyN9Yxug8uov7N4ancHPbhQxI+eW/Q9hHNytexqWcyICe7r8Youhd2ppNy+K0YxrkyC3nfKiEvNC7j/89Hmz+Vl2yYaTtZsCU6nV8Slar15IC0LIMTQ3QtKc7tjJ/zqVdvwJ79wUcP3UjHBtcam6TWrZbrTGkEvCeIVvxrAAUrpQUqpAuCLAG4N4H6boJRai5DdtPKKalzGp1wGQACzkFdRW3ps2Asm/Vj71YpywFIM+1DOutroCm0GlPwtcLB/4Npq7OalcPiBvWaNukciBM8/yyjy2bduMZzFU3vB8qP3PIv/+S8P4G/vebbtY9ZqFMeWy9gx7V7/AYxJTreM3c0jaDopeXZt2H/PzVY5KGIux8tOmPv+7OXW9+xaOwBrIXenLY9FWydWxJRkjiyXQKn/AqlfUrIARa81nbjs7rDtpLB+EkRg3wLA7vRz3Pxe4PzJ1560KtJssMZOvtp6EUUmLkLVactLfXvGPmiBvRSwFDNtSjAXbGzWIQ2jqdYDSn5OMvYP3KSHqyOD/W1hDylV1VrDYNUtl20CAERcevIlIQJVp5b97GK++cPsJQvaYXUGr8I+wIqRbhp7c8buXNrsRliWvYy4aLQ72q9iVNO90j7o9QJT7mKwid9OpU57ZxBgFDEPmoNQQSyfscPuz3nyYT7vWyfjbZPEftKz4ikh5HZCyF5CyN6FhYX2v+DCFnNS8NVXbsWhxSKOLpXwps8+aI0vt/LiAOqDCmdy3pdQDRn7ABZPg2p3BIAP/OIl+NhrLseVOyabfmZYw7p/8Go16run3p6xtx1Q6lHGrui1hl7yWdOUy+0zym7HioTs8tvuXvnUqVzbDzgbb9/eIrDHPNwSSy5F0IQUbVuLCFq6cxKXzGlZ2zEXXd4XO6aTOPzBV+HxP3kpXn3lVvzez50LAChUOxvZLzrqOrNp2ZpwDWJdpB0rsFecgd045m2TCah6eyO2fhFEYD8BYJvt31vN7zVAKb2DUrqHUrpndnZ2XQ90+w1n4/AHX4XLt01A0Wv4q2/twz1Pz+P/PnYSgPf2JMYmM1M4uepd9GjI2Adooa1eo1C0WmCTpwCwcyaJWy/f4jo9mpS9u2IKCvM8aR/YJ236fbatxt4rKUZvyNivPXsaV++cwltsQy0M50Jrdvn9d6+5ApduzeIXL9+Milpra/PKrjBbSjFCs2YNAGWluQiakKJQddqynS9fVQMPeHbcaiKtWnLTMREf+ZXLsGUiDkI6nxMpOeS/2ZRsSUHt6jedwuRcZ02Ovc5MThpUU8EgAvuDAHYTQnYRQiQArwHwtQDu1xOW9bCxXrZ+zWvfqfV75ofK7mPhhBkKERJ+gOkEL7/wsGi193StVN+/2Q57sbCdLGBJMQE/7x/9zrP4q7v3Wf+uqrUGv51MTMSXfvv5uGxbs0FZfaG1cUzMGGrzRAxfe9sLcNMFRi/7Qr51Ie2T9x8EAGydbJ2xu+08dTPz8uPn4za4EyRxl61PRktu6/doJEKQkoSOLbSLjmln++7YTtoa/cCSEKcR2NFl48rrgk0ZAINrA951YKeUagDeBuBbAJ4G8CVKabBrUxwwv2UmqbDugHybjH1zNoYNGbmhJctJvmK46CXE3gzL+IV9gOMhaqZ2ElLUsx2NvdnbZeBAe/nFjiXFBHyl9LHv7sfH733O+re93bEdzl2kLIgxDZl5q3uNlz9xYg1lRcdyUcHmbKxhhN1JTIy6DvyUVHcpBmg9RGZ3hAwDt9er1KKl0046JnRcPC05pp3fcN1O6+vZNgZznWLtPnVM9x5bLiElC9g5bUhxg+oWG0iUoJR+A8A3grgvPzh3RrInl/3fK0shhOCyrRP42Yk1z/vO2ZziBimwe61uC4uEJHhuUOoksDOZxzkE5UZPpRifzyPL2JnGzi79WWBnAcUtsOcqKm753/fj6l1TWCuruO3aHU23sSMLEdduL7d+dD/PVaHq3fobBExjb5BiFB2bJ9o/t17thK1wZuzbpxN47TXbQ5Gb6ivyGgN7oWrEh7qlwWBm7ENpKZCJCRAiBJrZqcDO/HUpxvvPOms2he8/M49ajbquaGOXr1qNWtrmIMBa26aS4WVgdpKS0fGgaLUmm2AW2P1m4w+896YmO1w3WKYXZMZuL2oaro6RDjN20yedZeyaDikasTJvVnh1k2KePml0cP30kCEZ7phONt3GjmfG7nJSj7nIIE4K1e4Gd9oRc83Y/XXipOTOMna9RlFRa033/YH/cYnv++gENiXtrJ0Y08OCdTIZ1Ix9KC0FCKkHdQDWJhWWAbTq1tgyGYeqUyx4XDqz8eREj3xL/HLILL5tn2odHIKCbVFyywg7ydgBw13Pz21Z4AryebcHj7WyCq1GQSk6lmJYxl5W9IZ1b5m4gJgYsTYG2TnoGDm/1HQ79CImRrFaVpt085Ji2DXbE5F2gV3Ta6ioNV/TwevFTee3u1C2IhUTrcXzfiiG7FTphC3ucA5RFczAzpLHTtxie8lQBnY7u2aS1uVrvqpZdqxebDHXoXmNA7OMvVeGVH5hi5DdhmjCwGvvaUXV8d2nzwDwH9j9EhNZ0Tq4D4t9OnOtrFqB0O+ykrrGbvxeVdMbrj4IIdiYiTW4FjLs35tIiDhnNtX2sRbyVVz7F99t+L5bQbJeuHTvimHae5iBMGXed8Em2RUVf66faVmwFo/4oRSyt7wbblcV7CqINWiMcldMX7lsa9bKIL3MgOxsmTC6Ek54tDwa9yEgKQkD5ceer6gQIsRyAQwbr72nd9x3EPc8PQ8geL2fEIK42L4/24kx1Vly/ZkzsNcXgvsLELLDdbGi1ppeg43ZmGvGvmS7KpxJya7Snx12wshVtAYjtLJLsGTH4PVcsUwzTCmGtQSyx9JM6c7PVUKnUkyvM3b2WM6MvWgG9qQURYTwjD1w3vOK85GUotg1k0JJMXYfrpbUtlkkG3I60TJjFwcuY89XNGTMZbu9wGvvqd3RLoxjia+jG+njP3gOL/zQ913nE+xdDSVFb9s55USKmgNKWr0rJuZYlOGVsS8Wqtg+lcCN58/h/be0t09qsC6wnZCKita0wNyt1dBO2IusgfqKPyaBljpoye20eNqPjD0pCQ1XI4BpvSwLIIQg3cKzvd8MbWD/7RedjSf+9GWWQ2GurDZscvciJQvIxkXXIKDpNZQUHemY4Guyr5fkzCuJXlHfotT44VPNTPKrv3N9KI8bX8fz/uFvPQMAeNyl28nu51FW9HrA8xkgmjP25unfTRNxnF6rNNkNLxUUbJ6I4VNvuAo3nNt+KM9ez7Bns85ResC9cGmn0IOMPWEeA9PKOwm+KdM91K9Fs5Wx92iOAzCO0XnFakgxxjEYi3t4xh44hBBrSGbNDOx+JtC2TMRdpRj7gBMrnuo1in+894BlW8D46qMn1rW3cb0YC3t70xED1C95ncXTZ07ncNP5c67DPEGQkJqtYP2y38VG156xV1S940xWirLiqXFMZZeMfctEHFqNYt6xXHqlpPhq82S88QW7rK+dZnTOLJgFdjdvGaB+pRVmxh6JECNAs8DuMiHrhZcXixdu1sVhk7T9bYDRYcW6YoDWW5b6zVAHdqBewFv1mbEDhhzjpsnavWbiooCyouM/Hj6OD939DP7sv+rLrZ84sYZ3fPFRvPlzewP6K9qTK/crY68HDk2v4bmFIs51MQ0LingLjxo37AZUz55pnii2OyCWVb3jTNbK2NW6xi47NHZmE/DcfGMXzGpJbZiObMelWyfwxduvBdA4bl9SdMQdUgzT2L2Kp60WZwdJUo7WpRiPvbFudBrY2RWT36J3ECTlaMPxVbUatBqtWwfHBKsjb9AY+sA+ba51Wy4oWCv788Y4ayaJ5xYKTT4bLLBnYgImEobR/mPHjF2N9n5o1sa297D3BKsXB+bzlkNlJ4Q9Hu7Ebe/pYkGBXqOhduakWnjUuGG/FHazilgpKvVCo6Jbl/R+B3fqGbtNY3cUjc+ZM7pdDi3WH59SitWyiskOJm8BuA6+lJXmzUvtpJheFE/Z/ReqjYHdz1UCK7z61dnZ8++3TTUIkg7raqfJYDYuYrU8mNa9Qx/YZ8zJv/l81SowtmPnTBI1iqZe9vrkqoidpnvf/QcM//eVYv2Si7U4rcfZ7Zf+8Ud4xcd+2PHv5Srhjoc7ccvYWedHq2URQTxuJxk764hKSFFXeW25qGCzeSIyMvbO2gBlS/IwXuuqVmsK7BvSMRDSOKSUrxqdLe3WATpJs4BnO7kVXaQYMRqBECGexVP7ezlM0rH6HuGiZXvhr3gK+Hd4ZEkYK2b3gqQsNNSYlq0hQeM1nUkbu08HkaEP7FPmpe6RZSOL9iPFMF9w57iw/YzM+t2Z3ap91NtuDNSJxvbkyTUrE6t1uNeTdev0CsuQyxZgWGDfEGJgT0rNLWatYCfZ8zamsVpSmwLdclHBxowReO0au28pxmECVlZ0xBxZI9Oa7Vl2aZ0ad33wpf6+8lqw0qo1tJ2FdVDYd4MyHdxvuyPgv11Q6YMUw+oH7Gqdta8ylYBS4/215DHs2E+GPrDHRCNzYcVNP2/kjNVJ0/imqi/qEJuClz2Y24cSfvvzD/k+1s//5Ij1tdfkqxt6jRoeFfHeSTFi1Jy8s2XPLDuZSXeWhXZC0qUToRXsdTnPdPdzZu0rJRVTSckIgmZXTIT478EXIgSE2NodNXdP/LSjL7uTQqKdlGOikVJq7kptvh+vjUuAkYjExWhoa/EY9gIiO2mGobEr5olVjPam3Rcw4kSN1o+R1WvYgpqrdhp7DJ5zWfrTb4Y+sBNCkIoJVh+xn0yMSRpOwyV7ljOXrgf2aIQ0BHN7kHcWzLxYK6u466Hj1rKP4yvuAzVuFKzj6l3GDjTvPWUdJu02IXX1mB0OrrDX4jpzS8/jxxtbHpcKVUyzwK4afexJSfDdg08IsfxlAHeN3Tpul+1bnXZxiNEIYmLEeg6qWg2UurcQxqWIqxRzJlfBv/zwUOCTwW6kYwKeWyhi7+Hljqyl2efUr8au6kbW3MuMfcK6sjfeYyxesOf1ws2GdW87y+Z+MPSBHTDeJKc7COzshXGOA9sDu/0NdN3Z08hXNUs+yVVUnDOXwk3nz/k2wnrqZA6qTvHmFxnLHDrZcJ5r41oZFs69pytFBQkp6svQa70kJcHaXesH9txcaPpj21sOVXPx8GRSsvZz2tvV/CJFI2aANYyonFIMYO7/tA0YdeOfn5LrWXCrTpOYx2KOHzxrbCj77Red1fFjdwq7Inj1P/3Y5hvvQ4rpNGPXey/FsASGJTROeYtp7csDuNR69AK7j+DnZbmZq6iQhIjl6HfF9glECHDljklQWp+sWyuryMQEbJmM48RKGf9y30FPiwLGIbOT5sXnGYMqXiPwbrDg1cviKdC893S5pISarbPHBJoHo7w4k6siQoyNNnEx2pA9sQ/kVFJCUo4iX1EN69cOWwBlMYqqVrOydjfL35gYaXBm7KT1z0nGplszecMt87cPc+UqKl704e/j/v2LOLZcQjRC8No2NsFBYF/1dyZX8W17wXT4faebZw/csNodQ5aW7NSte+uWJdEIsWQ89llot1S8H4xEYE/HBOuM7mddG8vq3TJ2+47Ou377Ojzy/pdayxTYh2y1pGIyIWHLRBz5qoYPfONpvP8/n2j5mCzgbJtMYCYld5Sx29swe4mzQ2W1pGIyZNtglk17ecE7ObVaxmxahixEMZOWGgO72ck0lZQwEZewVlZRqOpIdXiClIUIFK3WtGTDTlPGzjpExM5fs4aCpOpdkLRn7PvP5HFkqYQ/+b9P4uhyCZuysdD1dQD4f67bgZ+/bDMA4NFjq5hI+LO9YLbHdz10HPfvX2x7e0WrQYpGemapATRvUWItx+wYxGgE2bjYYP8wKIxEYLfLL34ydiEaMbsYmgO7XceORgiycbGp0LNaUpFNiFYbHVAftfdioVDBZEKEJESwbSqOYx1o7OwE1GuNPUKA+55dsOSNpWL4GTvLcP1eoi8VFevEm3Ss87Pa0xISMnERa2XVNHHqLIuWhAiqmm4VKt0Kr8a+0uaMfT0DQinbdqFWBcmYFEXZfEx2EjswX8CZXMXaMhY2shDF77zEkBcfPba6Ll3/0WPt50FUvdbTwingHdjtTKekBtuKQWEkArs94PkNfinHuDBQd3Z0wrIlVuhZNSUJZijGfrcVC/mqFYC2TiY6ytg79T8PClY8uvoDho3s6bVyqK2OQP0k7bczZqlQxbT5vDp9ZqzAnpIwa2bzzk33fpAFQ4qpZ+xuGnukYby/2MUIfNqmsZdbauz1x7RvcFoqdGZl0C32Pa6dTNp+6JcvBeCvq8Rt4UvYNAd21ZozYEwnJSxxjT0cds/Vfa79mgSlbJe7DK/pTiYPFKsaFK2GoqJjIi5iqy1jZ/3uXhgyhvGm3zoZx8nVsu8CYb8C+7tedp719UpRwXy+GrofPCu8+c3YFwsKZszn1ekzs1yqZ+yzKRlLRQX7Tuc7fh4lJsVobaQYtVmKWY8bod35sK7Vu3XF1E9k9sC+f74Q+pWVnZQsWDIh2yjlh1+9ahteuHum5XJ5Rj8Ce0yMQhIi1hVzziU+TCUlLsWExUVbMtbXfjU4Nz/o1ZLiWqC0SzH2tXCzaRn/38vPx80XzGG5qDRZFNix+9hsm0xA1amrh7cbubIKQnrfFXPBpgz+9Y1XAwCu+LPvgFJgm61YFgapDjR2SikWC1Vr+thp+btsXiJPJCRcsX3S+n6ni49lhxTjlbFXtGYpZj2e9XaNvT5b0Vpjd05ATvYwYwfqjQiXtNkS5WQ2JfsKjIre+8AOGMlUoxTTGB+mkv6Ov9d09UwRQj5MCNlHCPkZIeQrhJBwLP/acPGWzt5MQOOHh7FYUFw/9FbGrmjWtOpEQgIhBG958dm48fwNANDykszuzrjVlHD8dsbkKhrSstB2UUMYXL1rquHfZ8+Gu5ov0UFXTKGqoarVLNkhLgkNWfNioYps3KhrvOT8Oev7ncoUrHjKrgac7o7se40Zu7FfNbqO1ywdE1FQjPZaNkTnFtjjtisU5zLtqR5m7ABwuen2edMFc21u2Yi9ntAKRa/1pBjsJBsXG7pinA0M02bG3ukkedh0+0x9B8DFlNJLATwL4A+7P6TOYcNEnbSWpR1vqKqmY62sWjq4HVYA++tvP2u1Ldr712fMSbTFvPeZuyFjN7Nevzr7Wlm17Il7TUyM4sOvvtT6NxvKCAtrcMXHh511wLChr4QjYz+dqzT42tx+g9HX/QuXb+7omNiAEpNi3NsdjcDOxs/drHb9kpYFUGq0a/7x1540vudyJSnbrhIWC1VcZsuW/c5XBMXn33QNHvmjn8P5Gzt7f6TMwS67yZ4brCum1xSrGu5+8jS+8fgpz+JpjRrusoNEV9f2lNJv2/75EwCv7u5w1s+33nlDR9ppyjEpyCrbboGdGTkdXylbSx3s5k5MClgoVAA0Xz1oeg2FqmYd3+YJw7vEb2eMXzvisPiVPdvwwt2zULSa1eMfFnWPmvZSzLwV2I3g7SyenslVsMHWHfLeV16A977ygo6PSRaiULSaVah0k1fiUhQ1amSWshC1lh6vBxY8Hji0bH3PTYaIi8ZxaXoNC/kqztuYxmPm5G0vNXbAPCF3pnAZvxcToNWoq7maHeO91/vAziba73roOApVNymmPqTUy4J1O4J8pt4I4JsB3l9HnLcxjY0dtHjZp/uA+qUsy77t2D9U++eNQo89I2IFI6+M3epDN71eZCGKzdk4Hj66ird8/qG29gL9DuyAsddz+3S4+jrQWfHUytgztq4Yhxvlxsw6oo0DiWXsLTR25wkp18Vrxlp2nzxpBOlXXbLJ9XZM2stXNMznqg02GGHPGwRF2qcZWD+KpwDwtpecA8Bo/dVtXuwMp+3AoND2mSKE3EMIecLlv1ttt3kfAA3AnS3u53ZCyF5CyN6FhYVgjr4LUrHG1VwssHsV1v78Fy+2vo5GiBVM7L/jZezlNjl69lwK9z27gG8+cRov++h9rtt/GIMQ2HsFm+zz0+7IAju7yoqLUSi6kcGyLDaI9sx68dS7KybpWCXYzWvGssKfHV8DIcDf/NplrrdjFtVPnFxDvqrh3A1pqyus1xn7evFrLdAvjf2dN+/GZVuzVkumU15jfvsrwxbYKaU3U0ovdvnvqwBACHkDgFsAvJa2EMoopXdQSvdQSvfMzrbf/xg2rAjCPogs23aTYgDgtmt3WMFdr9EGSSImRpGWBU8zILfx8kttBd+iouONn33Q81gNC4PxCOwAMwJrL8WsOtpArQXPWg1LRQU1GozFMGt3bNXpwoq+7Dbd+OczKebxE2vYMhH3lL/Ye/iJE8bilrNnk/irV1+Ka8+awo7pcIvcQcEWbrSzalb71BUjRCM4d0Paqq055xKcfjKDQrddMS8H8G4Av0Ap9T9KOQCwAtiiGYwX2mTsAPCyizYCAC7e0lwgmk3Lnhm7m9OfU/Lxsl8FxitjB5o9arxg6wJZ5wl7fkuKZrWSznXY2ugGK55aa/VaDLGx19rY5rVOjd283F8tqdjZIkCzjP3Easn69y2XbsYXb3/+urpx+oFfX/Z+FU+BxgEsZ8bOmhqc9iT9pttn6u8BpAF8hxDyKCHknwI4pp7AFmk8edLIdhbyVaRkoWUBZzYt41Nv2IN/uu3Kpp/NpGTPjL3sMmTysos34vyNadxyqaGfrpVU14GliqpD0Wq+NkONCsZKMh+B3ZEVWxm7UsOZnPFaBJmxF6oaYmLEVRKwTirmcefKWtdSDACcu8F7vyz720+u+t9FMGj49WXvl8YOAOdtrA9AOgN7WjYSi0HL2LvtijknqAPpNVfvmgYA/PTQMq7ZNYXP/OiwLwMx1rPuZDYte+4ydVu6sCkbx93vvAEA8Pyzj+B9X3kCZ3KVBv8ZoH9Tp/0kKUcb7IK9yJUbVyGy4FpW9UC3PclCFFqNYrWkWNJB0zHbVgkqWg1lVV+3FGO/arx8u/doCLsiOGG2zfbaSygIklZ7a+uMt9rHwH7Ztvpr4Ez8CCGYsPW6DwojMXm6HlKygE3ZGD73kyO45i8NL5R8BwsenMykpKYhEQZrwfPyDdmcNYL5aZdJ1NxYBnbB14BSzjEwYpdi2Gsx7dLl1CmszW65qHhmxXWNXasvZFhnL7ldRrl8a6vAzqQYI7CHvbg6DPwu3FD1/kkxm7L1ZMutpTGb4IF9oGAfIFby/eAvXbLu+2JdNm7143be3KxNk3nK27E6asYpsPuVYspqY8Yu1jP25aKCbFwMpJOCZYpLRcUzeNo1dnaV1U3B+47XXYlXXLwR26a8vXlSkoAIMWSMpBQdGl3djrXj1UdXTL8ydqC+ku/s2VTTzybiIlbLIyTFDDvvfeUFeOudD1v/vrHDcWg7CUmA7jFoYQV2D29uNh15yjWw92Yp8SDhV4rJVxp1bEtjV3UsFRVMBzQwwrpSlouKdXXlhL0+a2U1kKusl160ES81i/VeRCIE6ZjhZTKMMgxgXA0JEdI2Y+9n8RQAvv72F6JQVV0ThcmE5Hq13U/GOmN/5SWb8KU3P2Yg9CEAABULSURBVN/6N9s+vh5SNgdIJ9bSBY+MfSIhQhYiOL3WbDFg7Tsdwsvs9WIs+PDfFcOoSzE6lgO0rmWZ4nJB8fT7T0hRRAjwsXv24+4nTwNAT5aPs8cY1hM/21k8yMVTwBiAvHLHlOvPuBQzgNgvdbu5lE1IjX3MdkqKDiFCPN+YhBBsysZw0iVjZ2/4Yc3I1kM2LpoGX95Ze61GUVAaR7wtKUbRMZ+veM4kdArT2PNVzfMESwhBjRoy0D//4KD1d4QNk3uGNbADzfYeTmo1Cq1G+zKg5IeJuDRwXTGD+Uz1kLl0DHNpGa++cmtX92N3gHRSUvS2CxemU7LlHGmH2R742Qw1KuyaSYJS4GgLj/uSqoPSxisZe1fM8ZVyS326E+wn5E78X+ZCXkoC2AP78J74U7LQUmPvxyLrTpjLyCgpuu8dAr1gMJ+pHhKNENz7rhc3OBiuB5axu2nDZR9Of15ZS6GigRDDuXBcOMu0Bm61gIE9V/YTHsvYT69VUNVqgW17sptPtZJX7nzTNXjlJXVdvBfyGTueXrs5Bkk61rpYzgJ7P0zA/MAkv5UB8mUfnzSwBevZcuOk1Uq3kqq3fYxUTHB1e8xXNaP7YQg7HtbLrhkW2L1XprG+Z3sGzQL7EdPnvtOFGl7YR/rtRltOrj9nBtefY2wE0mu0J4uXmdwzMcRdUylZwFKLoMgW2Axqxs6umpw7lPsJD+wBwQK3V/G03SaddIuMfZxkGMCQFebSMg4vegf2vEtRORIhkIWIJeG0CsKdYM+G/VgUuLXEhQU76Ux1UfjvN0lZaLlakgX2QdXY2VUTW4oyCAzmMzWEJOXWxVM/UozbScFrD+uosyET8/TeAeDp2xKXojiyZJwQ5gKw7AXQMA0c1H0GBXseNnVgWT1opGOtNXaVaewDGtidS68HgcF8poaQhG2k3Imf4qnTRphRqGpDOVHYLTMpyVp+4oalsTuem4QYtXr/gzAAAxq7W4K6CggKVkfo9Y7TIGnXFcOlmM4Zv4gRElbG7pJ5VFQdG9pkeilbV419YjFfXb+Z1DAznZKx77S3Rz3L8JyBPRUTgDVDbw/jhBiUbh8Ub7/xHGybiuOm89c/XNdvUrKIsqpD02sQXLLy+oKTwWwgYNPPg+TwOJinwCEkJkRBSIuMvc2b0sszI19Rx2o4iTGTkrFUUDx3YRar7hO5LIOdy8ihFC8HLbgkZQGvvWbHUBfXmZzmNW3Mds22+wz1i7QsgJD6lPggMH4RIyQiEWIsU3bJ2A0ppn1XDNBsX1qqrn8p8jAzk5Kg6DXD/tallY+dAJ195VtMPTwoGYbxnd+7YeDGxkcFaz1e1X1pO7O9dltJOAhEIgQpWRiojJ0H9gBJyIJrxl5R/RVPgeaFAyVFG9PAXl836BrYPbzRmRNfNh6s5rx7Qxq7W3ijc9ZPu/V45RYrCQeFTEwcKI19ME+BQ0pSat78QylFyU+7o8ebu6zqSIyhFMO6T+bz7llyvqq5eqNvMheotNjSyBkwWvksAWi5a3ZQyMRF3u44qiQkoUknLFQ11Gh7Lw8WpOwau6LVoOp0rKZOGSzzdrMyBoznye053Tpp/N5Ftp2ynMGGZexe6/EqbfYZDALZuDBQGfv4pYIh4rark23y2dimz7heQKr/ftllV+q40MrKGPBuA33+WdP43G9ejWvMDVmcwYdp7J5STIsl4oNCJibi6PLgrH3mgT1AEpKAVUcB5fSav92blsZue3OXVLZSb/xeprgUxZaJOJ44seb680LFPbATQvDC3bNhHx4nQJIeHWGMisbaHQdXYMjExbYLuXvJ4D5TQ4ixIMIjY28T2JNmVm5/c7fbvDTqXLljEo8cXXX9Wb6qdeS0yBlc2hZPWVeMMLifg0xMHKiuGB7YAyQhCU3tjqd9SjFCNIK4GG1Y6lse88D+vO0TOJ2r4LiLOVqhqo6l1cIowtYKemW8ZVWHLEQGulc/EzdsEZyT4/0ikMBOCPkDQgglhMwEcX/DSlKKNrU7nlgtYyIh+qroZ+JCg99EPWMfzwDGNtb87HizHOMlxXCGj2iEIClFPTP21ZIy8NPXzF3TbadCP+g6sBNCtgF4KYCj3R/OcJOQhabi6V17j2PbZMLX78+kZCza/FGKbVbqjTqsw8XZGVOr0bE1RxtVUjFvv5iVkhrYmsOwYDW0+by3cV0vCSJj/yiAdwMYjGuQPpKUolB1apkWLRaqUPSa77VZMykZC7Y3xrhLMRMJEZIQwRlHL/tKSYFWowPn28JZPynZe+/pakkZ+EUibFvWmQGZTu4qsBNCbgVwglL6mI/b3k4I2UsI2buwsNDNww4sTDJhWTvzmP6Dl57r6/enUxKWbQsHrEXWY5qZEkIwl5Yxn2vMgpid76A5LXLWTyomelr3rpRUTCYGPWM3koyhCeyEkHsIIU+4/HcrgPcCeL+fB6KU3kEp3UMp3TM7O5rtaEwyef9XnwQAq+h38WZ/wzKGFFO1pibZGz3tMmE5LrDnxA4L9DxjHx3SLbxWjIx9sAM7ey+eyQ2JFEMpvZlSerHzPwAHAewC8Bgh5DCArQAeJoRsbHV/owwL4F977CQAWHq53wA0nZRQ1WpWAdZtr+e44ZSngLqOGbTRF6d/bJ2Muw74UErNjH2wkxtZiGIqKQ1Mxr7uiEEpfRyAZQJtBvc9lNLFAI5rKLlkaxaXbs1aI9CrJQURggZ/9VawAtFyQTE1RxUJKYroALd5hc1sWsKjxxp72dnQUlDLqjn9Z8tEHMtFBRVVb+ggy1WMFsJBL54CRqIxNBk7pzMu3pK1MvXVkopsXPTdf1t3NDTO+uO6PcnOdFLGcrGKmtkfXFI0fOZHhwGMb7fQKMJM337nzocbvs/aBwddigHMHQLFEQvslNKd45ytM2ZTMlZKClSzG6aTos+2KaO979hyGYAxsDHugX0yKaFG62vH3vZ/HgEA3HLppn4eFidg2NXXd/fNN3x/yWwmmEoOthQDGDsEnPWgfjHeUSMEZtIyKAWWi4qRsXegDW6dTIAQ4LC5jLnIx+YxzeSpooKFfBXfMz/4H/mVy/p5WJyAufYsd9O2Y6buvtXnLEg/mTU7uCiloWzv6gQuxQTMLJNT8tWOM/aYGMXmbNxqk/SzBHvUYdrqUlHBA4eWAQCf+82rB9qbm9M5MTGKd73sPAB1m14AOL5iXL2yYbVBZi4dQ1WrDcSKPB7YA2Y6VQ9EqyW148GKHdMJK2Mv+9i8NOpsNhdn7D9TwKPHVpGQorj+7LF2rhhZMi6+7MeWS5hJSUNhq8HqBM4urn7AA3vA1KWDqtF/2+GKth3TyYaMfdwDO7sEf+9XHsddDx3HeRvTA20GxVk/abN7LG9bWHFitWztsR10WFuz19avXsIDe8BMJ+uDCkVF7zhj3zmdwHJRwVpZRVnRERcHP1MJk5gYhWTba+rXd4czfLjt/c1VNGQG3ACMwSahecY+gqRjAiKkXvTJdDhctGM6CQA4ulQa20XWThS9Zn192baJPh4JJ0zSLlJMoaL6ngPpN9ae3gHoZR/vdDAEIhGCbLy+JqvTbGNTtm4mxKUYg3+67Ur85OASful5W3zbM3CGj6TLirxhmuVIywJkITIQUsxwPGNDRjYuWtX8TrMN1kWzXFRQ1Wpj3xUDAC+/eCNefvHYOlWMDWz1nf0KLV/RhsZSgxCCuYw8ENa9XIoJgWxCsqSYTp0ZJ8xBjJNrxomBZ+yccUE2V99VzXZHvUZRUvShcjedS8cGQorhgT0EJuIiNHMEvlMpJi0LECIEJ1eNwD7Im9k5nCCRBSMcVc19BkySGRYpBjD8YgZBiuGBPQTsa7w6DeyEEEwkRJxcNd4c8SHo3+VwgoBl7IojsA9L8RQwAjvvihlRmJkXUO9r74SJhIQTq1yK4YwXstiYsbN+9mHR2AFjk1KuoqGi6qCU4l/uO4iDC4WeHwcP7CHApiUBrGv0fVM2hkOLxvTpuHvFcMYHNq9Q1Rz7CIboM2ANKeWqOLRYxAe+8TTe+JkHe34cw/OMDRG7ZpJd/f62qfoQTpJn7JwxIRIhEKOknrFXh2815OasMSV7cq1snZiWCv52HgfJ8DxjQ8QFmzJd/f5ZthMDz9g544QsRFFVmRQzfIF9i2lWdmKlbO0v7sThNSiG5xkbIjZlY9iYieHqXVPr+n17xp/kxVPOGCELERcpZniKp2zA8ORq2aqTmSuMewqPGiFACMF3fv+GdVvLMlsBAEjIXIrhjA9GYGddMUbxdJgy9pgYxUxKxpHlEu59xtgdwLZA9RJePA2JdEyEGF3f08s2KQHDVTjicLpFFqP1wF7RQMjwzXKcNZvEfc8uYLGgYCYloajo1mrHXsED+wDC+nmB9XXVcDjDiixEoJhSTL6qISUJQ2fTfPm2CctWgNXbSrblIb2Ap4MDyqd/4yrI68z4OZxhpUGKGSKfGDsXba43T5y3IY0f7l9EscdmZl1HDkLI2wkh+wghTxJCPhTEQXGAl5w3h+vO4ZuCOOOFLERRVozstqgMj7OjnRefN2d9fe6GNIBGx8pe0NWzRgh5CYBbAVxGKa0SQuba/Q6Hw+F4EZOiWCsbRdN8ZTiXuWfjIm65dBM0nVo7ews93oPa7bP2FgAfpJRWAYBSOt/9IXE4nHElJkQwb+rRhao2VB0xdv7+fz4PAPDj55YAAMUeZ+zdSjHnAnghIeQBQsgPCCFXBXFQHA5nPIlLUZTNwL5SVIZmLZ4XKZflIb2g7emQEHIPALctB+8zf38KwLUArgLwJULIWZQ2t+QTQm4HcDsAbN++vZtj5nA4I0pcjKKi6lD1Go6tlPGqSzf1+5C6ImnOoRSVAQvslNKbvX5GCHkLgC+bgfynhJAagBkACy73cweAOwBgz549fZjF4nA4g05MNIqn+88UoNeoVXwcVlhXT6819m6lmP8E8BIAIIScC0ACsNjtQXE4nPEkJkZR0WqWu+nQB3ZLihmuPvZPAfgUIeQJAAqA17vJMBwOh+OHmBiBotWszphh19jjYhQR0vviaVeBnVKqALgtoGPhcDhjDrMPYFuIhrUrhkEIQVISel485aONHA5nYGAWGvP5CggBUiPgbpqKCUPX7sjhcDiBEbcCe3UofWLcSMo8Y+dwOGNMTKoH9mGXYRg8sHM4nLEmZfZ9n1otIx0b7sIpIy0bUkxV03FkqQi9Bxa+PLBzOJyBYTppLoPOV5GJj0rGHkWxquNLDx7Diz58L+55+kzoj8kDO4fDGRhm0rL19ahk7EyK+aOvPgmgN735PLBzOJyBYSYlgZj10lHR2FM2jX3LRLxhp3FY8MDO4XAGBlmIYnPWWA05SoGdDVy9/rodPXlMHtg5HM5AsWWCBfbRkWIYWycTPXlMHtg5HM5AMZEwAvooZeyMaXPxRtjwwM7hcAaKuNnLnhnBjH2SB3YOhzOOJMzAPip+gvaMnV2NhA0P7BwOZ6C43lzivnvILXsZDYE93puMfTRELA6HMzK86pJNuOI9k1YRddhhW5QAQBJ6k0vzjJ3D4QwUhJCRCepAY8beK3hg53A4nBDJ9mFZCA/sHA6HEyKzNpuEXsE1dg6HwwkRQgj+8pcuwbYeDScBPLBzOBxO6Pz61dt7+nhciuFwOJwRo6vATgi5nBDyE0LIo4SQvYSQq4M6MA6Hw+Gsj24z9g8B+FNK6eUA3m/+m8PhcDh9pNvATgFkzK+zAE52eX8cDofD6ZJui6fvBPAtQshHYJwkruv+kDgcDofTDW0DOyHkHgAbXX70PgA3Afg9Sul/EEJ+FcAnAdzscT+3A7gdALZv722FmMPhcMYJ0o2DGiFkDcAEpZQSQgiANUpppt3v7dmzh+7du3fdj8vhcDjjCCHkIUrpnna361ZjPwngRebXNwLY3+X9cTgcDqdLus3YXwDgYzAknQqAt1JKH/LxewsAjqzzYWcALK7zd4cV/jePB/xvHg+6+Zt3UEpn292oq8DeDwghe/1ciowS/G8eD/jfPB704m/mk6ccDoczYvDAzuFwOCPGMAb2O/p9AH2A/83jAf+bx4PQ/+ah09g5HA6H05phzNg5HA6H04KhCuyEkJcTQp4hhBwghLyn38cTNoSQbYSQ7xNCniKEPEkIeUe/j6kXEEKihJBHCCFf7/ex9AJCyAQh5C5CyD5CyNOEkOf3+5jChhDye+Z7+glCyBcIIbF+H1PQEEI+RQiZJ4Q8YfveFCHkO4SQ/eb/J8N47KEJ7ISQKIB/APAKABcC+HVCyIX9ParQ0QD8AaX0QgDXAvidMfibAeAdAJ7u90H0kI8BuJtSej6AyzDifzshZAuA3wWwh1J6MYAogNf096hC4TMAXu743nsAfJdSuhvAd81/B87QBHYAVwM4QCk9SClVAHwRwK19PqZQoZSeopQ+bH6dh/GB39LfowoXQshWAK8C8Il+H0svIIRkAdwAw2cJlFKFUrra36PqCQKAOCFEAJDACDrDUkrvA7Ds+PatAD5rfv1ZAL8YxmMPU2DfAuCY7d/HMeJBzg4hZCeAKwA80N8jCZ2/BfBuALV+H0iP2AVgAcCnTfnpE4SQZL8PKkwopScAfATAUQCnYHhMfbu/R9UzNlBKT5lfnwawIYwHGabAPrYQQlIA/gPAOymluX4fT1gQQm4BMO/HlmKEEAA8D8DHKaVXACgipMvzQcHUlW+FcVLbDCBJCLmtv0fVe6jRkhhKW+IwBfYTALbZ/r3V/N5IQwgRYQT1OymlX+738YTM9QB+gRByGIbUdiMh5PP9PaTQOQ7gOKWUXYndBSPQjzI3AzhEKV2glKoAvozx2eVwhhCyCQDM/8+H8SDDFNgfBLCbELKLECLBKLZ8rc/HFCqmFfInATxNKf2bfh9P2FBK/5BSupVSuhPG6/s9SulIZ3KU0tMAjhFCzjO/dROAp/p4SL3gKIBrCSEJ8z1+E0a8YGzjawBeb379egBfDeNBut2g1DMopRoh5G0AvgWjiv4pSumTfT6ssLkewOsAPE4IedT83nsppd/o4zFxguftAO40E5aDAH6jz8cTKpTSBwghdwF4GEbn1yMYwQlUQsgXALwYwAwh5DiAPwbwQQBfIoT8JgyH218N5bH55CmHw+GMFsMkxXA4HA7HBzywczgczojBAzuHw+GMGDywczgczojBAzuHw+GMGDywczgczojBAzuHw+GMGDywczgczojx/wOGHN8G75jBVAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "N_sim = 10000\n",
    "e = np.random.randn(N_sim)\n",
    "te = np.arange(N_sim) * Ts\n",
    "\n",
    "_, d, _ = control.forced_response(Hu, te,e)\n",
    "plt.plot(te,d)\n",
    "np.std(d)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAD8CAYAAABpcuN4AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl8XGXd9/HPb2ayNVvbLF3SNk13UpZCYxcWWQpSFCi3oLSsKogLPKg8LuV5bkVRHkXvWxRvUCsFqyJtBZSICLKv3dICpTvpmqalzb7v+T1/zGkZYjKZNJmcyczv/Xr11Zkz17nmd1jyzXWuc64jqooxxhjTE4/bBRhjjIlsFhTGGGOCsqAwxhgTlAWFMcaYoCwojDHGBGVBYYwxJigLCmOMMUFZUBhjjAnKgsIYY0xQPrcLGAiZmZk6ceJEt8swxpghZePGjeWqmtVbu6gIiokTJ1JUVOR2GcYYM6SIyP5Q2tmpJ2OMMUFZUBhjjAnKgsIYY0xQFhTGGGOCsqAwxhgTlAWFMcaYoCwojDHGBBUV91EMJR/UNPPqrqMcrW0hIyWBs6dkMiFjmNtlGWNMjywoBkl9Szs/+ed2/rzuAJ1dHlN+4UnZ/J9PnsSkrBR3ijPGmCAsKAbB0bpmPvfwBnZ8UMv183K5bl4uEzKGcai6mb+/e4jfvb6HS3/1Bj/+9CksmpXjdrnGGPMRFhRh1tTawRdXFLG3vIHln/sY50/PPv5ZXmYyty+YymcLxnP7Y2/ztZXvUF7fyk1n57lYsTHGfJRNZofZ9wu3srm0hvuXnP6RkAg0Oj2RP948h0tOHs0Pn97GH9eGtPyKMcYMCguKMHr9/TJWFZXwpY9P5qL8UUHbJvi8/GrJ6SyYkc1dT23h5R1HB6lKY4wJzoIiTNo7Ornrqa1Mykzm6xdODWkfn9fD/UtOJ39sGrc/9jYllY1hrtIYY3oXUlCIyEIR2SkixSKytJvPE0RklfP5OhGZGPDZnc72nSJycW99isjrIvKO8+eQiPytf4fojic3lbKnvIGll8wgMc4b8n7JCT5+fe1sELjtsbdpbe8MY5XGGNO7XoNCRLzAA8AlQD6wRETyuzS7CahS1SnAfcC9zr75wGJgJrAQeFBEvMH6VNVzVHWWqs4C1gBP9v8wB1d7Rye/fPF9Ths/vNdTTt0ZP3IY9155Ku+WVPPLF3eFoUJjjAldKCOKOUCxqu5R1VZgJbCoS5tFwArn9ePAAhERZ/tKVW1R1b1AsdNfr32KSBpwATDkRhTPbztCaXUTt543Gf8/hr775CljuPKMcfz21T1sO1Q7wBUaY0zoQgmKHKAk4P1BZ1u3bVS1HagBMoLsG0qfVwAvquqQ+yn5yFv7GDciiQUn9X00Eeg/P3USw4fFsfTJzXR0vUvPGGMGSSRPZi8BHuvpQxG5RUSKRKSorKxsEMsKbucHdazfW8kN83Pxek5sNHHMiOR47rpsJpsP1vCHNfsGpD5jjOmrUIKiFBgf8H6cs63bNiLiA9KBiiD7Bu1TRDLxn576R09FqeoyVS1Q1YKsrF6fDT5o/vp2KV6PcOUZ4wakv0tPHcM5UzO57/ldVDW0DkifxhjTF6EExQZgqojkiUg8/snpwi5tCoEbnddXAS+pqjrbFztXReUBU4H1IfR5FfC0qjaf6IG5obNT+fu7h/j41EwyUhIGpE8R4buX5tPQ2sF9L9jEtjFm8PUaFM6cw23Ac8B2YLWqbhWRu0XkcqfZciBDRIqBO4Clzr5bgdXANuBZ4FZV7eipz4CvXUyQ006Rqmh/FaXVTQO+XtO0UalcO3cCj647wK4jdQPatzHG9Eb8v/gPbQUFBVpUVOR2GXzvqS2sLiph439eRHLCwC6jVdnQynk/e5nZuSN45PNzBrRvY0xsEpGNqlrQW7tInsweUlSVF7cf5ZypWQMeEgAjk+P58nmTeXlnGRv3Vw54/8YY0xMLigGy60g9pdVNLJjR/cJ/A+FzZ04kMyWB/3rO5iqMMYPHgmKAvLjjCADnhzEohsX7uPX8yazZU8GbxeVh+x5jjAlkQTFAXtp+lJNz0hiVlhjW71kyZwJj0hP52XM7iYb5JWNM5LOgGAC1zW1sOlDFedPCN5o4JjHOy+0LpvJOSTWv7IqcGw2NMdHLgmIAFO2rpFPhzMkZg/J9V54xjrHpiTz4cvGgfJ8xJrZZUAyAdXsqifd6OH3CiEH5vnifh1s+PokN+6pYv9eugDLGhJcFxQBYu6eCWeOHkxQf+nMn+uvqj00gIzmeB2xUYYwJMwuKfqprbuO90hrmTRo5qN+bFO/lC2fn8equMraU1gzqdxtjYosFRT8V7auiU2HepMGZnwh0/fxcUhN8PPiKjSqMMeFjQdFPRfsr8Xpk0OYnAqUlxnHDmbn8c8sHFB+tH/TvN8bEBguKfnq3pIYZo1MHdX4i0OfPyiPO6+HhN/e68v3GmOhnQdEPnZ3KuwermTV+uGs1ZKYk8OnTc3hi40Eq7XkVxpgwsKDohz3lDdQ1t3Oai0EBcNPZebS0d/KntftdrcMYE50sKPrh3ZJqAE53OSimjkrlvOlZ/GHNPprbOlytxRgTfSwo+uGdkmpSEnxMykpxuxRuPnsS5fWtFL57yO1SjDFRxoKiH949WM2p49LxesTtUjhrSgYzRqey/PW9tligMWZAWVCcoNb2TrYfruWUcelulwL4n6198zmT2HmkjtfftyXIjTEDJ6SgEJGFIrJTRIpFZGk3nyeIyCrn83UiMjHgszud7TtF5OLe+hS/e0Rkl4hsF5Hb+3eI4bGnvJ62DiV/TJrbpRx32WljyEpN4KE37FJZY8zA6TUoRMQLPABcAuQDS0Qkv0uzm4AqVZ0C3Afc6+ybDywGZgILgQdFxNtLn58DxgMzVPUkYGW/jjBMth+uBYiooEjweblxfi6v7Sqj+Gid2+UYY6JEKCOKOUCxqu5R1Vb8P7gXdWmzCFjhvH4cWCAi4mxfqaotqroXKHb6C9bnV4C7VbUTQFWPnvjhhc/2w3XE+zzkZSa7XcpHLJkzgXifhxVv2aWyxpiBEUpQ5AAlAe8POtu6baOq7UANkBFk32B9TgauFpEiEfmniEwN7VAG1/bDtUwblYLPG1nTPBkpCVx26lie2HSQ2uY2t8sxxkSByPop55cANKtqAfA74OHuGonILU6YFJWVDf6T3rYfrmPG6Mg57RToc2dOpLG1g8eLDrpdijEmCoQSFKX45wyOGeds67aNiPiAdKAiyL7B+jwIPOm8/itwandFqeoyVS1Q1YKsrKwQDmPglNW1UF7fwkkRND8R6JRx6ZwxYTh/WLOPzk67VNYY0z+hBMUGYKqI5IlIPP7J6cIubQqBG53XVwEvqf9i/kJgsXNVVB4wFVjfS59/A853Xp8L7DqxQwufHR/4J7JPGp3qciU9u/HMieyraOTV9+252saY/uk1KJw5h9uA54DtwGpV3Soid4vI5U6z5UCGiBQDdwBLnX23AquBbcCzwK2q2tFTn05fPwGuFJH3gB8DNw/MoQ6cY1c8zYjQEQXAJSf7L5Vd8dY+t0sxxgxxvlAaqeozwDNdtn0v4HUz8Jke9r0HuCeUPp3t1cCnQqnLLe8fqSczJYGRyfFul9KjeJ+Ha+dO4BcvvM/e8oaIuzrLGDN0ROJkdsTbXVbPlOzI/8F7zdwJxHmFP6zZ53YpxpghzIKij1SV3WUNEbEQYG+yUxP55CljeLzoIA0t7W6XY4wZoiwo+qiyoZWapjYmD4GgAP+kdl1LO09usktljTEnxoKij3aXNQAwOSvyTz2B/1kZp45LZ8Wa/baqrDHmhFhQ9NHusnqAITOiEBFunD+R4qP1vFlc4XY5xpghyIKij3YfrSfB5yFneJLbpYTs0tPGkJEcz4o1+9wuxRgzBFlQ9NHusnomZaXgiYCHFYUqwedlyZwJvLj9CCWVjW6XY4wZYiwo+mh3WcOQmZ8IdO28CYgIf1pnq8oaY/rGgqIPmts6KKlqHDLzE4HGpCdx8cxRrNpQQnNbh9vlGGOGEAuKPthf0YgqTBqCIwqAG+ZPpLqxjcJ3DrldijFmCLGg6IP9Ff5LYydmDM2gmJs3khmjU/n9W/vsUlljTMgsKPrggDMRnJsxzOVKToyIcMP8iWw7XMvG/VVul2OMGSIsKPrgQGUjqYk+0pPi3C7lhF1x+ljSEn383laVNcaEyIKiD/ZXNJKbMQz/48CHpmHxPj5bMJ5nt3zAkdpmt8sxxgwBFhR9UFLZyISRQ/O0U6Dr5+fSocqj6w64XYoxZgiwoAhRR6dSUtXIhJFDcyI7UG5GMudPz+bP6w7Q2t7pdjnGmAhnQRGiD2qbaevQqBhRgH9V2fL6Fv655bDbpRhjIpwFRYiOXRo7VK946uqcKZnkZSbbpLYxplchBYWILBSRnSJSLCJLu/k8QURWOZ+vE5GJAZ/d6WzfKSIX99aniPxeRPaKyDvOn1n9O8SBcWyNpGgZUXg8wg3zc3n7QDWbD1a7XY4xJoL1GhQi4gUeAC4B8oElIpLfpdlNQJWqTgHuA+519s0HFgMzgYXAgyLiDaHPb6nqLOfPO/06wgGyv6IRn0cYk57odikD5srZ4xgW72XFW7b+kzGmZ6GMKOYAxaq6R1VbgZXAoi5tFgErnNePAwvEfw3pImClqrao6l6g2OkvlD4jyoHKRnJGJOHzRs/ZurTEOK48Yxx/33yIivoWt8sxxkSoUH7q5QAlAe8POtu6baOq7UANkBFk3976vEdENovIfSKSEEKNYRctl8Z2dcP8XFrbO1m5oaT3xsaYmBSJvx7fCcwAPgaMBL7TXSMRuUVEikSkqKysLOxFlVY3MW7E0HlYUaimjkrlrCkZPLp2P+0ddqmsMebfhRIUpcD4gPfjnG3dthERH5AOVATZt8c+VfWw+rUAj+A/TfVvVHWZqhaoakFWVlYIh3Himts6KK9vZWx69AUFwI3zJ3KoppkXth9xuxRjTAQKJSg2AFNFJE9E4vFPThd2aVMI3Oi8vgp4Sf3LkxYCi52rovKAqcD6YH2KyBjnbwGuALb05wAHwuEa/1IXY4fQ40/7YsFJo8gZnmSXyhpjutVrUDhzDrcBzwHbgdWqulVE7haRy51my4EMESkG7gCWOvtuBVYD24BngVtVtaOnPp2+HhWR94D3gEzgRwNzqCfuUHUTADlReOoJwOsRrp+fy9o9lez8oM7tcowxEcYXSiNVfQZ4psu27wW8bgY+08O+9wD3hNKns/2CUGoaTKXHgiJKRxQAVxeM577nd7FizT7+33+c4nY5xpgIEomT2RHnUHUTIjAqLXruoehqRHI8i2aN5a+bSqlpanO7HGNMBLGgCMGh6iayUxOI90X3P64b5k+kqa2DvxTZpbLGmA9F90++AXKoujlqJ7IDnZyTTkHuCP64dj+dnfaoVGOMnwVFCA5VN8VEUADccOZE9lc08uqu8N+bYowZGiwoeqGqlFY3RfVEdqCFM0eTnZrAijX73C7FGBMhLCh6UdnQSkt7J2OjaDHAYOJ9Hq6ZO4FXdpaxt7zB7XKMMRHAgqIXh6qj+2a77lwzdwJxXmGF3YBnjMGColfH7qGIpaDITk3kU6eM4S9FJXaprDHGgqI3h2LgZrvu3HzOJBpaO1i14YDbpRhjXGZB0YvDNU0kxnkYPizO7VIG1ck56cybNJJH3txHm60qa0xMs6DoxZHaFkanJeJfozC2fPGcSRyuaeaZ9w67XYoxxkUWFL04UttMdhQv3RHM+dOzmZSVzO9e34N/MWBjTCyyoOjFkdrmqF7jKRiPR7jp7Dy2lNaybm+l2+UYY1xiQRGEqnKktoVRqRHxNFZXXHnGOEYMi+Oh1/e6XYoxxiUWFEHUtbTT1NYRsyMKgMQ4L9fPy+XFHUfYU1bvdjnGGBdYUARxtNZ/s112WuyOKACunz+ROK+H5W/YqMKYWGRBEcSR2hYgup9DEYqs1AT+Y1YOT2w6SGVDq9vlGGMGmQVFEB84z8qO9aAAuOmcPJrbOnl07X63SzHGDLKQgkJEForIThEpFpGl3XyeICKrnM/XicjEgM/udLbvFJGL+9Dn/SLi6knxI3XHgiK2Tz0BTBuVyrnTslixZh/NbR1ul2OMGUS9BoWIeIEHgEuAfGCJiOR3aXYTUKWqU4D7gHudffOBxcBMYCHwoIh4e+tTRAqAEf08tn47WttCaqKPYfEhPVo86n353MmU17faE/CMiTGhjCjmAMWqukdVW4GVwKIubRYBK5zXjwMLxH8r8yJgpaq2qOpeoNjpr8c+nRD5GfDt/h1a/8XyPRTdmTdpJKdPGM5vX9tDuy3rYUzMCCUocoDAXyEPOtu6baOq7UANkBFk32B93gYUqmrQdSNE5BYRKRKRorKy8DyNzR8UdtrpGBHhq+dN4WBVE3/ffMjtcowxgySiJrNFZCzwGeBXvbVV1WWqWqCqBVlZWWGpx3+znY0oAi2Ykc20USn8+pXd9lxtY2JEKEFRCowPeD/O2dZtGxHxAelARZB9e9p+OjAFKBaRfcAwESkO8VgGVGencrQudtd56onHI3zlvMnsOlLPSzuOul2OMWYQhBIUG4CpIpInIvH4J6cLu7QpBG50Xl8FvKT+VeQKgcXOVVF5wFRgfU99quo/VHW0qk5U1YlAozNBPuiqGltp61A79dSNy04dy7gRSTz4SrEtFmhMDOg1KJw5h9uA54DtwGpV3Soid4vI5U6z5UCG89v/HcBSZ9+twGpgG/AscKuqdvTU58AeWv/YzXY983k9fOnjk9h0oNoWCzQmBoR03aeqPgM802Xb9wJeN+OfW+hu33uAe0Lps5s2KaHUFw7H7qHIjuEFAYP5TMF4fvni+zz4ym7mTcpwuxxjTBhF1GR2JCmv848osiwoupUY5+ULZ+fx2q4ytpTWuF2OMSaMLCh6UF7vX9MoM8WCoifXzcslNdHHr1563+1SjDFhZEHRg/L6FpLivCQn2F3ZPUlLjOMLZ+Xx3NYjbD1kowpjopUFRQ/K61vITI13u4yI94Wz80hN9HH/izaqMCZaWVD0oKK+lYxkO+3Um/QkG1UYE+0sKHpQXt9i8xMhslGFMdHNgqIH5fUtZNmpp5DYqMKY6GZB0Y2OTqWyodVGFH1gowpjopcFRTcqG1rpVLs0ti9sVGFM9LKg6EZ5vf9mOwuKvjk2qvjFCzaqMCaaWFB0o+L4zXY2R9EX6Ulx3Hz2JJ7fdoS3D1S5XY4xZoBYUHTj2Igiw0YUfXbTOXlkJMdz77M7bGVZY6KEBUU3jgVFlgVFn6Uk+Ljtgims3VPJ6++Xu12OMWYAWFB0o6y+hXivh7QkW77jRFwzdwI5w5P46XM77Cl4xkQBC4pulNe1kpESj4i4XcqQlODzcsdF09hSWsszW4I++twYMwRYUHTD7sruvytOz2H6qFT++1+7aOvodLscY0w/WFB0o6Khxa546ievR/jWxdPZW97AX4oOul2OMaYfLCi64T/1ZCOK/lpwUjazc0fwixd20dja7nY5xpgTFFJQiMhCEdkpIsUisrSbzxNEZJXz+ToRmRjw2Z3O9p0icnFvfYrIchF5V0Q2i8jjIjKoj0NVVWdEYUHRXyLC//nkDI7WtfCbV/e4XY4x5gT1GhQi4gUeAC4B8oElIpLfpdlNQJWqTgHuA+519s0HFgMzgYXAgyLi7aXPb6jqaap6KnAAuK2fx9gntU3ttHWonXoaILNzR3LpqWNY9tpuDtc0uV2OMeYEhDKimAMUq+oeVW0FVgKLurRZBKxwXj8OLBD/JUOLgJWq2qKqe4Fip78e+1TVWgBn/yRgUK+vrGz035U9MtmCYqAsvWQGnQo/fXan26UYY05AKEGRA5QEvD/obOu2jaq2AzVARpB9g/YpIo8AHwAzgF91V5SI3CIiRSJSVFZWFsJhhKaywR8UIywoBsy4EcO4+ew8/vp2Ke+UVLtdjjGmjyJyMltVPw+MBbYDV/fQZpmqFqhqQVZW1oB9d5UTFCOHWVAMpK+eP4XMlAR++PQ2W9rDmCEmlKAoBcYHvB/nbOu2jYj4gHSgIsi+vfapqh34T0ldGUKNA8ZOPYVHSoKPb108jY37q3h6s92EZ8xQEkpQbACmikieiMTjn5wu7NKmELjReX0V8JL6f20sBBY7V0XlAVOB9T31KX5T4PgcxeXAjv4dYt9U2amnsLlq9njyx6Txk3/uoKm1w+1yjDEh6jUonDmH24Dn8J8KWq2qW0XkbhG53Gm2HMgQkWLgDmCps+9WYDWwDXgWuFVVO3rqExBghYi8B7wHjAHuHrCjDUFlYyvxPg/J8d7B/NqY4PUId12WT2l1Ew+8XOx2OcaYEIW06p2qPgM802Xb9wJeNwOf6WHfe4B7QuyzEzgrlJrCpaqhlZHDbJ2ncJk7KYNPn57Dstf28B9n5DA5a1BvkzHGnICInMx2U2VDm512CrM7P3kSCXEe7npqq01sGzMEWFB0UdXYysjkOLfLiGpZqQl88xPTeaO4nH+8ZxPbxkQ6C4ouqhpaGWGXxobddfNymTk2jR8+vY36FlsHyphIZkHRRWVjq10aOwi8HuGHV5zMkdoWfvH8LrfLMcYEYUERoL2jk5qmNhtRDJIzJoxgyZzxPPzmXjYftDu2jYlUFhQBapraULWb7QbT0ktOIjMlgW8/vtkecGRMhLKgCFDVaDfbDbb0pDh+dMXJ7Pigjt++utvtcowx3bCgCFDZ0AbYOk+D7RMzR/OpU8dw/4vFFB+tc7scY0wXFhQBPlw51i6PHWw/uHwmwxK8fPvxzXR02r0VxkQSC4oAVbYgoGsyUxK467J8Nh2oZsVb+9wuxxgTwIIiwPERhZ16csUVs3K4YEY29z67g+Kj9W6XY4xxWFAEqGpoZVi8l8Q4WxDQDSLCTz59CsPivXxj1Tt2FZQxEcKCIkBlo92V7bbstER+/OlTeK+0hvtffN/tcowxWFB8RFWD3ZUdCRaePIarZo/jgZeL2bi/0u1yjIl5FhQBKhtt5dhIcddl+YwdnsQ3Vr1Lg60FZYyrLCgC+J9FYZfGRoLUxDh+/tlZlFQ18r2ntrpdjjExzYIiQHVjK8NtjiJizMkbyf+6YCpPbDrIX4pK3C7HmJhlQeHo6FTqWtpJT7IRRST52oKpzJs0ku8+tYVdR+yubWPcEFJQiMhCEdkpIsUisrSbzxNEZJXz+ToRmRjw2Z3O9p0icnFvfYrIo872LSLysIgMyk/uumb/goDD7dRTRPF6hPsXn05KQhxffXQTja02X2HMYOs1KETECzwAXALkA0tEJL9Ls5uAKlWdAtwH3Ovsmw8sBmYCC4EHRcTbS5+PAjOAU4Ak4OZ+HWGIqhv96zzZiCLyZKcl8svFs9hdVs93/2bzFcYMtlBGFHOAYlXdo6qtwEpgUZc2i4AVzuvHgQUiIs72laraoqp7gWKnvx77VNVn1AGsB8b17xBDU9PkDwobUUSms6ZkcrszX/HndQfcLseYmBJKUOQAgTOJB51t3bZR1XagBsgIsm+vfTqnnK4Hnu2uKBG5RUSKRKSorKwshMMIrrrJRhSR7vYFUzl3WhZ3FW6haJ/dX2HMYInkyewHgddU9fXuPlTVZapaoKoFWVlZ/f6yamdBwPQku+opUh2brxg3Yhhf/tMmDlU3uV2SMTEhlKAoBcYHvB/nbOu2jYj4gHSgIsi+QfsUkbuALOCOUA5iINTaqachIX1YHL+7YTbNbR186Y8baW7rcLskY6JeKEGxAZgqInkiEo9/crqwS5tC4Ebn9VXAS84cQyGw2LkqKg+Yin/eocc+ReRm4GJgiaoO2qpwNpk9dEzJTuW+q2fxXmkNdz75Hv7/1Iwx4dJrUDhzDrcBzwHbgdWqulVE7haRy51my4EMESnGPwpY6uy7FVgNbMM/13Crqnb01KfT12+AUcAaEXlHRL43QMcaVE1TG8nxXuK8kXw2zhxzUf4o/vdF0/jr26Xc/2Kx2+UYE9V8oTRS1WeAZ7ps+17A62bgMz3sew9wTyh9OttDqmmgVTe12WhiiLntginsrWjgvhd2kTMiiatmD8oFcsbEHFd+KEei6sY20m35jiHF//yKUzlS28zSJzYzOi2Rs6dmul2WMVHHzrM4apvaGG4jiiEn3ufh19fNZnJWCl/+00a2H651uyRjoo4FhaO6qdVOPQ1RaYlxPPL5j5Gc4OVzj6ynpLLR7ZKMiSoWFI7qxja7NHYIGzs8iRVfmENLeyfXPLSWD2qa3S7JmKhhQeGoscnsIW/G6DRWfH4OVQ1tXPvQWirqW9wuyZioYEEBNLd10NLeSbqNKIa808YPZ/mNBRysauKGh9cfX8PLGHPiLCgIWBDQlu+ICnMnZfDb62ez60gdNyxfR02jhYUx/WFBgd2VHY3Om57Ng9fOZvvhOq55aC2VDa1ul2TMkGVBwYcLAtpkdnS5KH8Uy26YTfHReq753VrK6mzOwpgTYUHBh6eebEQRfc6bns3Dn/sY+yoaWLxsjV0NZcwJsKDAnkUR7c6aksmKz8/hg5pmrvz1WxQftWdvG9MXFhRwfLLTTj1Fr7mTMlj1pfm0tHdy5a/X2IOPjOkDCwr8p568HiElwZa+imYn56Tz5FfOZGRyPNc+tI5nt3zgdknGDAkWFHy4fIf/Md8mmk3IGMYTXzmTk8ak8ZVHN/LQ63vseRbG9MKCAqhparf5iRgyMjmex744j4vzR/Ojf2znm3/ZbE/KMyYICwr8p57SLChiSlK8lwevPYOvXziVJzYdZPGytRyptSuijOmOBQX+JcbTEm1+ItZ4PMLXL5zGb647g11H6rj8f96wSW5juhFSUIjIQhHZKSLFIrK0m88TRGSV8/k6EZkY8NmdzvadInJxb32KyG3ONhWRQXkKTV1zG2mJNqKIVQtPHsMTXzmTBJ+Xq5et5cFXiunstHkLY47pNShExAs8AFwC5ANLRCS/S7ObgCpVnQLcB9zr7JsPLAZmAguBB0XE20ufbwIXAvv7eWwhq21uJy3JRhSx7KQxaTx9+9ksnDmanz67k8//foOtPmuMI5QRxRygWFX3qGorsBJY1KXNImCF8/pxYIH4LyFaBKxU1RZV3Qs2TfBMAAANH0lEQVQUO/312Keqvq2q+/p5XH1S19xGqo0oYl5aYhz/c83p/OiKk1mzp4JP3v86r79f5nZZxrgulKDIAUoC3h90tnXbRlXbgRogI8i+ofQ5KFrbO2lu67Q5CgP4n8N93bxc/vbVs0hJ8HH98vX837++R0NLu9ulGeOaITuZLSK3iEiRiBSVlZ34b311zf67sm1EYQLlj03jH7efw81n5/Hn9QdY+MvXWLunwu2yjHFFKEFRCowPeD/O2dZtGxHxAelARZB9Q+kzKFVdpqoFqlqQlZXVl10/orbZ/5uizVGYrhLjvPznpfmsumU+grDkd2v57t+22MOQTMwJJSg2AFNFJE9E4vFPThd2aVMI3Oi8vgp4Sf23uxYCi52rovKAqcD6EPscFMdHFAk2ojDdm5M3kn9+7RxunD+RR9ftZ8F/v8rf3i61O7pNzOg1KJw5h9uA54DtwGpV3Soid4vI5U6z5UCGiBQDdwBLnX23AquBbcCzwK2q2tFTnwAicruIHMQ/ytgsIg8N3OH+u9qmYyMKCwrTs+QEH9+/fCZP3Xo2OcMT+fqqd7j2oXW2Eq2JCRINvxUVFBRoUVHRCe37zHuH+eqjm3j26+cwY3TaAFdmolFHp/Ln9Qf46bM7aGztYMmc8Xz9wmlkpiS4XZoxfSIiG1W1oLd2Q3Yye6DYZLbpK69HuH5eLq988zyumzuBx9aXcN7PXuGBl4tparU1o0z0ifmgOH7qyS6PNX2UkZLADxadzL++8XHOnJzBz57byXn/9TIr3tpniwyaqBLzQVHX3IYIJMdbUJgTMzkrhWU3FLDqlnnkjkzmrsKtnPuzl3nkzb0WGCYqxHxQ1Da3k5rgw+OxZ1GY/vE/RW8ef/7iXHIzkvnB37fx8Z++zG9f3W2X1JohLeZ/ja615TvMABIRzpycyZmTM1mzu4JfvfQ+P/7nDn754vt8tmA8XzgrjwkZw9wu05g+saBoardLY01YzJ+cwfzJGWw9VMPyN/by6Lr9rFizj0/kj+LaubmcPSXTRrJmSIj5oPAvCBjz/xhMGM0cm87PPzuL7yycwR/W7OPP6w7w3NYj5AxP4uqPjeezBeMZnZ7odpnG9Cjmf0LWNreTMzzJ7TJMDBiVlsi3Lp7B7Qum8q+tR1i54QA/f34Xv3hhF+dNz2bRrLFceNIokhNi/n9LE2Fi/r9I/0OLUt0uw8SQBJ+Xy04by2WnjWV/RQOrNpTw5KZSXtpxlKQ4LwtOyuby08Zy7vQsEnxet8s1xoKi1p6XbVyUm5HMtxfO4JufmM6GfZUUvnuIZ947zNObD5Oa6OP86dlcmD+Kc6dlkW7/nRqXxHRQdHYq9S3tNkdhXOfxCHMnZTB3Ugbfv3wmbxSX84/Nh3l5x1EK3z2EzyPMnTSSBTNGcd70LPIyk/E/G8yY8Ivpn5ANre10Kva8bBNR4rwezp+ezfnTs+noVN4pqeKF7Ud5YdsR7n56G3c/DWPTEzlzSiZnTcngrMmZZKfZZLgJn5gOijrnWRQ2ojCRyusRZueOZHbuSL6zcAb7Kxp4/f1y3tpdzgvbj/D4xoMATMlO4WMTRzI7dwRnTBhuIw4zoGL6J2StsyCgzVGYoSI3I5ncjGSum5dLZ6ey7XAtbxSXs2Z3BU9vPsRj6w8AMDI5njMmDOf0CSM4dVw6M8emMzI53uXqzVAV00FhIwozlHk8wsk56Zyck86Xz51MZ6dSXFbPxv1VbNpfxcYD/lNWx4xJTyR/TBozx6aRPzad/DFpjBuRZDf9mV7F9E/IWmf9HZujMNHA4xGmjUpl2qhUlsyZAEB1YytbD9Wy9VAN2w7VsvVQLS/vPEqn8xiaxDgPkzJTmJydwpSsFCZnJzMlO4WJGckkxtmlucYvpoPCRhQm2g0fFs9ZUzI5a0rm8W3NbR3s+KCOHYdrKT5aT3FZPe+UVPH05kMce46ZCIxJS2TcyGGMHzGM8SOTnL/9r0elJtpIJIbE9E9Im6MwsSgxzsus8cOZNX74R7Y3tXawp7ye3WUN7D5aT0llIyVVjbxZXM6RumYCH4YZ5xWyUxMZlZbA6PREslMTGZ3ufz8qLZFRaYlkpyaQkuCzSfUoEFJQiMhC4JeAF3hIVX/S5fME4A/AbKACuFpV9zmf3QncBHQAt6vqc8H6FJE8YCWQAWwErlfV1v4dZvdsRGHMh5Livcwc65/47qqlvYPSqiZKqpo4UNlIaVUTR2qbOVLbzM4P6nhtVzn1Le3/tl+818OI5DhGJicw8tjfwz76fviwOFITfaQmxpHm/B3vi/knIESUXn9CiogXeAC4CDgIbBCRQlXdFtDsJqBKVaeIyGLgXuBqEckHFgMzgbHACyIyzdmnpz7vBe5T1ZUi8hun718PxMF2VdvURoLPY8skGNOLBJ+XSVkpTMpK6bFNfUv78fA4WtvC0bpmKhvaqGpopaKhlcqGFrZU11BR30Jt87+Hyke/z0Na0kcDJC0xjpQEH0nxXoY5f5LifR++jvMyLN5HUryHpLiA7fFeEnxe4rxio5sTFMqv0nOAYlXdAyAiK4FFQGBQLAK+77x+HPgf8f8bWQSsVNUWYK+IFDv90V2fIrIduAC4xmmzwuk3LEFRY8t3GDNgUhJ8pGSlMDlImBzT1tFJVWMrlQ2t1DS2UdfcTm2z/++65jZqA/6ubfJvL61uor65nabWDhrbOujo1F6/p6t4n4cEr4d4n/9Pgu/D1/HHt3uPbz/W1ucVfB4PPo/g9Yr/b4+HuID3Po+/ndcjxHk8eD1y/H3gvnEeDx4PeEXweASPCB4Bj/jbivjvn/Fv9392/L3HeS/+0PN6hPSkOLxhni8KJShygJKA9weBuT21UdV2EanBf+ooB1jbZd8c53V3fWYA1ara3k37AVfV2MqIYRYUxgy2OK+H7FT/3MaJam3vdEKjncbWDppaO2hq63Be+7cd297a0UlLeyet7Z20tHfQ6rxu7Ti2rfP4tpqmtn9r19GptHX4/253/pxIUIXDC3ecy5Ts3sO5P4bsyXkRuQW4BWDChAkn1Mep44aTlxnef8DGmPA4NhJIx51f9lT1o8HRobR3dn7kfVunEy4d/raB71WVDlU61b/uXKfzuqOz+8/823G2q7MdslISwn6soQRFKTA+4P04Z1t3bQ6KiA9Ixz+pHWzf7rZXAMNFxOeMKrr7LgBUdRmwDKCgoOCEov3W86ecyG7GGIOI/9RSLExxhnJpwQZgqojkiUg8/snpwi5tCoEbnddXAS+pqjrbF4tIgnM101RgfU99Ovu87PSB0+dTJ354xhhj+qvXEYUz53Ab8Bz+S1kfVtWtInI3UKSqhcBy4I/OZHUl/h/8OO1W45/4bgduVdUOgO76dL7yO8BKEfkR8LbTtzHGGJeIamRMyPRHQUGBFhUVuV2GMcYMKSKyUVULemtnd7UYY4wJyoLCGGNMUBYUxhhjgrKgMMYYE5QFhTHGmKCi4qonESkD9p/g7plA+QCWMxTYMccGO+bo19/jzVXVrN4aRUVQ9IeIFIVyeVg0sWOODXbM0W+wjtdOPRljjAnKgsIYY0xQFhTOwoIxxo45NtgxR79BOd6Yn6MwxhgTnI0ojDHGBBXTQSEiC0Vkp4gUi8hSt+sJNxF5WESOisgWt2sZDCIyXkReFpFtIrJVRL7mdk3hJiKJIrJeRN51jvkHbtc0WETEKyJvi8jTbtcyGERkn4i8JyLviEhYV0WN2VNPIuIFdgEX4X/k6gZgiapuC7rjECYiHwfqgT+o6slu1xNuIjIGGKOqm0QkFdgIXBHl/44FSFbVehGJA94Avqaqa3vZdcgTkTuAAiBNVS91u55wE5F9QIGqhv2+kVgeUcwBilV1j6q2AiuBRS7XFFaq+hr+54XEBFU9rKqbnNd1wHbC+Az2SKB+9c7bOOdP1P82KCLjgE8BD7ldSzSK5aDIAUoC3h8kyn+IxDIRmQicDqxzt5Lwc07BvAMcBZ5X1ag/ZuAXwLeBTrcLGUQK/EtENorILeH8olgOChMjRCQFeAL4uqrWul1PuKlqh6rOwv/M+TkiEtWnGUXkUuCoqm50u5ZBdraqngFcAtzqnFoOi1gOilJgfMD7cc42E0Wc8/RPAI+q6pNu1zOYVLUa/zPoF7pdS5idBVzunLNfCVwgIn9yt6TwU9VS5++jwF/xn04Pi1gOig3AVBHJE5F4/M/5LnS5JjOAnInd5cB2Vf252/UMBhHJEpHhzusk/Bdr7HC3qvBS1TtVdZyqTsT///FLqnqdy2WFlYgkOxdoICLJwCeAsF3NGLNBoartwG3Ac/gnOVer6lZ3qwovEXkMWANMF5GDInKT2zWF2VnA9fh/w3zH+fNJt4sKszHAyyKyGf8vQ8+rakxcLhpjRgFviMi7wHrgH6r6bLi+LGYvjzXGGBOamB1RGGOMCY0FhTHGmKAsKIwxxgRlQWGMMSYoCwpjjDFBWVAYY4wJyoLCGGNMUBYUxhhjgvr/OZ4/OGDrygwAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Define a discrete-time power spectrum with sample time Ts = 1e-3, bandwidth 10 rad/sec, and \n",
    "# standard deviation 3 (when filtered by white noise)\n",
    "\n",
    "wu = 2 # bandwidth of the force disturbance\n",
    "std_du = 3\n",
    "Ts = 1e-3\n",
    "Hu = control.TransferFunction([1], [1/wu, 1])\n",
    "Hu = Hu * Hu\n",
    "Hud = control.matlab.c2d(Hu, Ts )\n",
    "t_imp = np.arange(5000)*Ts\n",
    "t, y = control.impulse_response(Hud,t_imp)\n",
    "y = y[0]\n",
    "plt.plot(t,y)\n",
    "std_tmp = np.sqrt(np.sum(y**2))#np.sqrt(trapz(y**2,t))\n",
    "Hu = Hu/(std_tmp)*std_du\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3.391379569767212"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl83GW9//3XNZN9su/71ibd0j3dKGUrVTYpCqgoKuA54HFD8b651bOIR8/veNSj56gocCuiggoiKCKI7G2htE33dMvW7NtkJvuemev3xyShpWmbZJbvLJ/n49GH7WT6/X4S6jtXru91fS6ltUYIIUTwMBldgBBCCM+SYBdCiCAjwS6EEEFGgl0IIYKMBLsQQgQZCXYhhAgyHgl2pVSiUupppdRJpdQJpdQmT1xXCCHE3IV56Dr/C/xNa32LUioCiPHQdYUQQsyRcneDklIqATgEFGvZ7SSEEIbzxIi9CLACv1RKrQT2A/dqrQfPfJNS6m7gbgCLxbJ28eLFHri1EEKEjv3793dprdMu9j5PjNjLgXeAzVrrPUqp/wX6tNb/er6/U15erisqKty6rxBChBql1H6tdfnF3ueJh6fNQLPWes/kn58G1njgukIIIebB7WDXWrcDTUqpRZMvbQWOu3tdIYQQ8+OpVTFfAJ6YXBFTB9zpoesKIYSYI48Eu9b6EHDReR8hhBDeJztPhRAiyEiwCyFEkJFgF0KIICPB7gWHm3r408EWnE7ZiCuE8D1PrYoRwITDyff/XsVDb9YC0DM0xh2biwyuSggRamTE7kFP7GnkoTdr+XB5LityE3hkRx3SPkcI4WsS7B4yMu7gZ2/Usr4wmf+6eQWf3FRIa+8Ih5t7jS5NCBFiJNjdNOFworXmF7tO0943wpeuLkEpxbYlGYSZFC9WthldohAixMgcuxvsg2Nsf3AXTfZhAK4ty+SShakAJMSEc3lpGk9XNPOlraVER5iNLFUIEUIk2OdgcHSCH71azeiEk+Otfeytt09/rCAlhm/fVHbW+z9zxQJufWg3j71dzz9dscDX5QohQpQE+xx88y/HeKqiGYDU2Ai2lKTywdU5XFaaRkSYifio8LPev64wmUsXpvL4Ow185vJilFJGlC2ECDES7LP04tE2nqpo5vNXLuQfthQRFW4mKvzi0yvXLc/i688epdY6yML0WB9UKoQIdfLwdBZ6hsb4//54hJW5Cdx7dQmJMRGzCnWASxakAPBOnc2bJQohxDQJ9ll47nArfSMT/McHlxNuntuXrCAlhsz4KPactl/8zUII4QES7LPw9P5mlmbFU5aTMOe/q5RiQ3Ey79TZZLOSEMInJNgv4lR7P0eae7llbe68r7GxOAVr/yh1XYMXf7MQQrhJgv0iXjrWjlLwgZXZ877GxmKZZxdC+I4E+0XsrLZSlp1AWlzkvK9RmBJDUkw4R6W9gBDCByTYL6BvZJwDjT1cVprq1nWUUizLTqCyVYJdCOF9EuwXsLvWhsOpuawkze1rLcuJp6p9gLEJpwcqE0KI85Ngv4AdVVYsEWZW5ye5fa2y7ATGHE6qO/s9UJkQQpyfBPsF7KzuYtOCVCLC3P8yLcuOB+BYS5/b1xJCiAuRYD+P+q5BGu1Dbs+vTylMsRAVbqKqQ0bsQgjvkmA/j5ePdwBwean78+sAJpOiMMUia9mFEF7nsWBXSpmVUgeVUs976ppGOd01yH+/fIr1hckUpFg8dt0FabHUWQc8dj0hhJiJJ0fs9wInPHg9Q4w7nPzrnyoxK8VPPrbao9cuTrPQ1D0sK2OEEF7lkWBXSuUC1wM/98T1jPTr3Q3squniX29YSnp8lEevXZxmweHUNNplOkYI4T2eGrH/D3A/cN6hqFLqbqVUhVKqwmq1eui2nqW15ok9DazJT+Sj6/M9fv0Faa5+7DWdEuxCCO9xO9iVUjcAnVrr/Rd6n9b6Ea11uda6PC3NMw8kPa3BNkSddZAPrpl/w68LWZgei0nBiTZZ8iiE8B5PjNg3AzcqpeqB3wNXKaUe98B1fW7vZM/0jUXJXrl+TEQYC9JiqWyR1gJCCO9xO9i11l/TWudqrQuBjwKvaa1vd7syA+ytt5MUE+7VI+yW5yZwRIJdCOFFso59ktOp2VltZUNRilcPnV6ek4C1f5TOvpE5/93+kXFu+dnb/MOv9uFwyqEdQoiZeTTYtdZvaK1v8OQ1feVgUzcdfaNcU5bp1ftMncJ0dB6j9h1VXVQ0dPPKiU5ePdHh6dKEEEFCRuyTXjzaToTZxFVL0r16n6VZ8Sg1v2CvaLATblbERobxwtE2L1QnhAgGYUYX4A+cTs2Lle1sKUklPircq/eyRIZRnGqZ1wPU/Q3drMlPIjcphpePtzPucM75cG0hRPCTVACeP9pGS88w21fn+OR+izPjqbXObS370NgEx1r7KC9MYtvSDPpGJnj4zVovVSiECGQhH+yjEw6+99JJlmTFc/3yLJ/csyjVQqN9iHHH7FsLHGrqweHUlBcks21pBjesyOL7f6/i3/9ynJFxhxerFUIEmpCfivnzoVaa7MM8dmcZZpP3VsOcqSh1qrXA0PRu1IvZX98NwJr8JMwmxQ8/sorU2Egefes0IxMO/s8Hl3uzZCFEAAn5EfvTFc0Up1k81p53NgpTXR0j6+fQwvdISy/FqRYSYlzPAMLNJh64cRmf2FjAU/uaaO+d+/JJIURwCulgb7IPsbfezs1rcr26dv29iieD/fQcgv1EWx9LJk9hOtPdlxXj1JrfvFPvqfIuSGstUz9C+LmQDvYXK11LBm9cme3T+yZZIkiMCZ91sPeNjNPcPczSrHODPS85hm1LM3hiT6NPAvfbfz3B2m+9LG0RhPBjIR3su2ttFKdZyEuO8fm9C1Mssw72k22u4/SWZMXN+PFPbiqkZ2icl461e6y+mQyNTfDbPY0Mjjn44ctVXr2XEGL+QjbYHU7NvvpuNhSlGHL/4lTLrOfYp7pBLplhxA6wqTiFvORonqpo8lh9M3nhaDvD4w5W5yfyZpUV++CYV+8nhJifkA32OusAA6MTrC1IMuT+RakWWntHGB67+PTJibY+EmPCyTzPwR8mk+LWtXm8VWOjyT7k6VKnPbmvkaJUC9++qYwJp+Yvh1u9di8hxPyFbLAfbnbNEa/MTTDk/tMrY2wXH7WfaOubbEVw/ge8t6zNRSn4g5dG7bXWAfbVd/ORdXksy05gSVY8fzzQjNbSjEwIfxO6wd7UgyXCTPEs15F7WnGaK9hrOi98uPWEw8nJ9v7zTsNMyU6M5rKSNH6/r4ne4XGP1TnlqX1NmE2KD61x7c792IZ8jjT38natzeP3EkK4J2SDfX9DN6snN/sYYWF6LGaT4mT7hU9TqrcNMjrhnHFFzHt9eVsp9sExvvX8cU+VCbgO+P7jgWa2Lk4nPc41HfTh8lzS4yL56Rs1Hr2XEMJ9IRns/SPjnGzvM2x+HSAyzMyCNMv0ipfzOdDQA7gO6LiYVXmJ3HFJIc8ebKG523Nz7buqu+gaGOPD5XnTr0WGmfn0pUW8VWPjeKsc9SeEPwnJYD/Y2INTQ3mhccEOrt7sh5p6LjhPvbOmi/S4SEpmearTXZcWobXmwdc9N5J+s8pKZJiJS0tSz3r9I+vyiAwz8fieBo/dSwjhvpAM9oqGbkwKVucbG+wbi1OwDY5R1THzPLvWmrdqurh0Yeqsd8ZmJ0ZzxyVF/G5vk8c2Ee2strK+KJmocPNZryfGRHDjymz+dLCF/pGz5/VPtPVx35OH+PqzR6m1Xvg5ghDCs0Iy2A80dLMkK57YSGN7oG0qdq2h313bNePH621D2AfHWD/Hw7Xv3VpChNnEMwda3K6xzjpArXWQrYtnPoDk9o0FDI05ePbgu/fSWvP//OEwfzrUwm/3NPLZxw/I6hkhfCjkgl1rzbHWXlYYtMzxTHnJMeQlR7O7buaVJQcbXR0d5/qTRUJMOFcsSuP5I61un4368nHXEXxXL82Y8eMr8xJZnpPAE+80Tof3zuoujrX28Z8fWs5/37qSUx397Kye+ZuXEMLzQi7YuwbG6B4apyR95u35vrauMJn9Dd0zjmgPNvYQGxnGwlnOr59p+6ocOvtHefs8Pw3M1svHO1iaFU9u0vnbLny4PJdTHf3TU0pP728m2RLBTatzuGFlFmlxkfxi12m36hBCzF7IBXt1h2sVSmmGfwR7eUEyXQNjNNjOXcVysKmblXkJ81qSuXVJOvFRYTzrxnRM18Ao+xu72Xae0fqUa8qyMCn465FWRsYdvHqig/cvyyAyzExkmJnbNxTwZpV1Tt0shRDzF3LBXjUV7JnGbEx6r6mVOfvq7We9Pjzm4ERbP6vz5veANyrczA0rs3mxsp2+kfltWHrtRCdac9FgT4uLZGNxCs8fbeONU1YGxxxcW/buaVS3rc8jzKT4zW5ZPSOEL4ResHcOkBgTTlpspNGlALAwLZbU2EjeOGU96/U9p22uo/DcWJL5kfI8hscdPHfowj1d7INj7D1tP+f1vx/vICcxmmUz9IF/r+tXZFFnHeQ/XjhOiiWCSxa821wtPT6K61dk8du9DV7tZSOEcAm9YG/vpzQ9zqcHa1yIyaS4tiyTV092MDg6Mf36qyc6iQ43s7F4/t0nV+S6err8bm/jBd9331OH+PDDu89639DYBDurrWxbmjGrr9U1yzIBaLIP89H1eYSZz/6ndf81izEpxTf/cuys11t6hr3SAkGIUOZ2sCul8pRSryuljiuljiml7vVEYd6gtaaqo5+SDP+Yhplyw4osRsad090Sazr7eXp/M1ctTj9n7fhcKKW4bX0ex1r7ONo885r2ffX26Z8WHtlRN/0Q92+V7YxOOLlulgd8p8RG8tkrFnD98iw+e8XCcz6ekxjNvVtLeOVEJ3+rbOMvh1tZ+62X2fyd19j632/QOMMzBiHE/HhixD4BfEVrvRTYCHxOKbXUA9f1uM7+UfpGJvzmwemUdYXJrM5P5IG/HOMHfz/Fp39VQWS4iW98wP0v4/ZVOUSEmfjjgeYZP/7MgRYsEWa+dVMZp7sG2VXjWkXzh4pm8pNjWDeHqaD7r1nMgx9fg+U8+wPuurSIpVnxfPaJA3zhdweJCjdzz+XFjI47+cZzlXP/5IQQM3I72LXWbVrrA5O/7wdOADnuXtcbqvxsRcwUk0nx0O1rKS9I5kev1WBWip99fC3p5+m/PhcJ0eFctSid54+0MeFwnvUxp1PzyokOLl+Uxi1rcslLjuZrzxzlb5Xt7K6z8ZF1eR6dsgo3m3j4E2u5fkU2Vy1O58+f38zXrl3CP15WzOunrDL/LoSHeHSOXSlVCKwG9szwsbuVUhVKqQqr1freD/vEqfapYPevqRiAjPgofvPp9fz5c5t56cuXsWmB50522r4qm66B0XM2Qh1p6cXaP8q2pRlER5j50UdX0z04xmce30+KJYJPXVLosRqm5CXH8OPbVvPoHetInXyAPTXds6PamH8XQgQbjwW7UioW+CPwJa31Oe3+tNaPaK3LtdblaWlpnrrtnFR3DJBiiSDFT1bEvJdSipV5iYSbPftM+8rF6SRbInhkR91Zr798vB2zSXHlIle7gNX5STx5zyauWJTG925d4bOWCwvSLGQnRElvdyE8xCMJopQKxxXqT2itn/HENb2hqtP/Hpz6QlS4mc9esYCd1V3sPiM8XzneybrCJBJjIqZfK8tJ4LE713PV4guvXfckpRRrCpI42NDts3sKEcw8sSpGAb8ATmitf+B+Sd6htaamc8BvWgn42u0bC8iIj+Te3x+kvmuQRtsQpzr62bY00+jSAFiTn0Rr7wgdfSNGlyJEwPPEiH0z8AngKqXUoclf13nguh5lGxyjf2Ri+ki6UBMVbuY7N69geMzBJx/dO33y0bYlvhuZX0hZjqspmxzaIYT7PLEqZpfWWmmtV2itV03+esETxXnSVJ+SqUOkQ9GVi9J57K719AyN8ft9TawvSiY/5fzNvXxpSZbrJ6njbRLsQrjL2IbkPnTa6gr24hAOdoC1BUm8cO8Wfr7zNJ++tMjocqbFRYVTkBLDsVbPHA4iRCgLmWCv6xokzKTISYw2uhTD5SbF8MCNy4wu4xxLs+JlKkYIDwiZXjGN9kHykmPO6WEi/Mey7HjqbUPnHLMnhJibkEm5+q4hCvxkPlnMbHluIgBHPXRWqxChKiSCXWtNo32IgmQJdn+2YnJlzJHzNCwTQsxOSAS7fXCMgdEJ8lNC+8Gpv0uyRFCQEsPhph6jSxHC42wDo/zLn476pE11SAR7w2RzqUKZivF7K3ITZcQugo7Dqbn1od08ua+JivpzD7XxtNAIdptrqaPMsfu/5TnxtPQMYx8cM7oUITxmV00XdV2DfPeWFWz1wabAEAn2IZRyLfMT/q0s2zXPLuvZRbBwODU/fLmKFEvErA+ucVdIBHujbYis+Ci3TiMSvrFsMtgrW2Q9uwgOv3zrNIeaevi3DywlMsw3GRQSwd5gH/KbrfPiwhJiwslLjqZSRuwiCNR3DfK9l05x9ZJ0blyZ7bP7hkaw2wYpSJYVMYFiaVY8J6RnjAgC//niCSLMJr5903KPnkZ2MUEf7AOjE3QNjFGQKiP2QFGcFkuTfQiHUxtdihDzdry1j5eOdXDXpUVkJrh/zOVcBH2vmEaba6mjjNgDR2FKDOMOTWvPMHmyqYxd1V0kRIezPDfB6FLELNgGRvnxazW8caqT+Kgw7jKg2V7wB7tdljoGmoLJjWT1tsGQD/Zd1V3c/gvXEcJfvrqUL25d6NMf6cXsjU44ePD1Wp450Exz9zAA/3XzchKiw31eS9AHe8PkiF0engaOwulgH2JLicHFGOzB12tIiA5nXWESP3ylihcr2/j6dUu4rNSYc4PF+T1zoIUfvVoNwH/fupLLStNIizPmfOWgn2Ovtw2RbIkgPsr33zXF/GTERxIVbqJh8nCUUFVrHWB3nY27LyvmodvX8v1bV2LtH+VLTx5iZNxhdHniDFprfvV2PSXpsey8/0puXptrWKhDCAR7g22Q/BD/cT7QKKUoTLFQbwvtYH/+cBtKwa1rcwkzm7hlbS4//fga7INjPLmvyejyxBl219o42d7PP2wp8ovpw6AP9urOAUrSY40uQ8xRQUoM9ZPTaKHq5RPtrMlPIj3+3RUV64uSKS9I4sev1dDeKwd/+4tHdtaRFhfJ9lU5RpcCBHmwdw+OYe0fpTQjzuhSxBwVplhotIXukkfbwCiVLX1cuejsuXSlFN+6qYyB0XG++9JJg6oTZ+rsG+HNKiu3rc/3m93tQR3sVR39AJRkyIg90BSkWBhzOGnvC81R6Vu1NgAuLTn3IemSrHhuWZvL80faGJ2QuXaj/fVoG1rDDSt80wdmNgI22LXWVHX003eBY9Smgn1RpozYA03h5Iay+hB9gLqr2upau54z89r1zQtSGZtwyhmxfuC3expZnZ/oVzMDARnsEw4ntz60m/f9cAdb/ut1ajoHZnzfqY5+4iLDyIz37a4v4b7CM9ayhxqtNbuqu7hkQQpm08xr1lfnJwFwSA4lMVRn/wjVnQNcV+Y/o3UI0GB/7nArFQ3d/OOWIpSCT/9qH+MO5znve7vWxqr8RNnQEYAy46OICDNN70MIJXVdg7T2jnBpSep535MRH0l8VBh11tD7xudP9td3A7C2MMngSs7mkWBXSl2jlDqllKpRSn3VE9c8n3GHkwdfr2FxZhxfv24J/3XzChpsQ7xwtO2s99VaB6izDrJtqfeb2gvPM5kUBckxnA7BqZhd1V0AbFl4/k1ISikKU2VJqNH21XcTGWaaPkfAX7gd7EopM/AgcC2wFLhNKbXU3evOZGhsgnt+s59a6yD3bStFKcW2JRksSLPw8Jt1aP3uCopXT3QA+OS0EuEdBSmW6dOvQsnO6i7yk2Muulta1vobb3+DnZV5iUSE+dfkhyeqWQ/UaK3rtNZjwO+B7R647jm+/sxRXj/VyTdvXMb7lmUCrpHdPZct4HhbHzsnRzpaa5450EJZTjw5idHeKEX4QFFqDA22IZwhtORx3OHknTrbBadhphSmWmjpHpaVMQYZGpvgWGsf5QX+NQ0Dngn2HODMbXDNk6+dRSl1t1KqQilVYbVa53Wjz19VwpN3b+JTlxSe9fr21dlkxEfy0Ju1gOvHo5Pt/dy+oWBe9xH+oSDFwuiEk47+0FnyeLiph4HRCbYsnEWwp8Tg1Ew3nBK+daiphwmnptzP5tfBhw9PtdaPaK3LtdblaWnza2C0MD2W9UXJ57weGWbmk5sKebvWRoNtkP95pYqkmHC/2QUm5mdqZUwozbM/e7CFiDATlyy4eLBPdcEMxekqf7D3tB2lYG3+uZlkNE8EewuQd8afcydf86mbVrtC/NaHdvN2rY2vvG8R0RH+sQtMzM9Uq+X6rtBYGdMzNMYfDzTzwVU5JMRcvGldUerkktAQ+fr4m921NpZlx8/qv5WveSLY9wElSqkipVQE8FHgOQ9cd05yEqP54lULGRl3cM/lxXx8Q76vSxAelp0YTbhZ0WAPjRHpwzvqGBl3cuelhbN6f1JMOHFRYTJiN8DIuIODjT1sKk4xupQZud2PXWs9oZT6PPASYAYe1Vofc7uyebjvfYu4732LjLi18AKzSZGXHDN9ClYwGh5z8POddTR3D/PU/iZuXZvL4sz4Wf3dd7tgBu/Xx1/tb+hmzOGc1ZSZETxy0IbW+gXgBU9cS4gzBXtwffelk/zyrXqiwk3csiaXb25fNqe/X5ASQ2VLr5eqE+ezu9aG2aRYN8MzP38Q9CcoicCWnxzDnjobWuug20Fs7R/l17sb+NiGfP7jprJ5fX6FKRZerGxn3OEk3Oxfa6mD2e46G8tzEoiN9M8IlX8Jwq8VpsQwOOaga2DM6FI87qVj7TicmjsuKZz3N62ClBgcTk2LLHn0mYHRCQ439bBpgX/Or4MEu/BzBanBu6Rvd52N7IQotw6CKUwN3WZpRtlRZWXCqbnCj8+dlWAXfm1Bqiv0aq0zd/AMVFpr9tTZ2Fic4tYU09SS0FBslmaUl493kBQTzlo/3HE6RYJd+LWcpGiiwk1UdwRXsNd0DtA1MMZGN5fLpcVGEhNhlhG7jwyMTvDKiQ6uWpxBmB8/0/DfyoTAteRxQVosVefpuR+o3p48IcndYFdKTTZLkxG7L/x6dz39IxN86hL/blciwS78Xkl6LDWTp2EFi7dqushLjr5oB8fZKEyJkRG7Dzicml+/3cCWklRW5CYaXc4FSbALv1ecFktr7wjDY8HRxdDp1Oytt7OxyDOrKgpTLTTZh5iY4bAZ4Tk7qqy0943wsfX+v6tdgl34vWBb+VHdOUDP0PiMDe3moyjVwrhD09IjSx696fkjbcRFhQXEGQ8S7MLvFU2dfxokXR731tsBPBrs4DpST3jHuMPJKyc62LYkw+8O1ZiJ/1coQl5hqmse+nSQjNj3nbaTHhdJfrL78+vwbrCflvNPvWbvaTu9w+O8vyzT6FJmRYJd+L24qHBSYyOCZsReUW9nXVGyx1okpFgiiIsKC6m+9b72t8p2osJNXFbiv5uSziTBLgJCYYolKPqOt/QM09o7wjoPbm5RSlGcapFg9xKnU/PSsXauKE0PmDMeJNhFQChMtQTFVMy+0675dU93BSySYPea3XU2OvtHuSZApmFAgl0EiKJUC9b+UQZGJ4wuxS376u3ERYbNuuf6bBWlxtLSM8zQWGB/ffzNvno7n33iADmJ0Wxb6v+rYaZIsIuA8O4xcIE9Kt1Xb2dNQRJmk2dbEJdmuHrq1ATZDl0jdfaNcPevK0iKCeeXd67D4qctemciwS4CwtTB1oG8lr13aJyqjgHKvdA8qiQjDoCqIOupY6Q/Hmihe2icRz5ZTunk1zdQSLCLgDC15DGQR+xHJ086Wp3v+WAvTIkhwmyiOshaLxjpxco2VuYmBFyogwS7CBAxEWFkxEdyOoBXxhxp6QFgeU6Cx68dZjZRnGahSoLdIzr7RzjS3BtQ8+pnkmAXAcN1/mkAj9ibeylIiSEhJtwr1y/NiJOpGA/ZWdUFwBWL0g2uZH4k2EXAKEq1BPRJSkdbeinzwmh9SmmGa2VMoK8c8gc7qq2kxkawNMuzq5d8RYJdBIyCFAtdA2P0j4wbXcqc9Q6N09w9TFm294J96gGqzLO7x+nU7KzuYktJGiYPr17yFQl2ETCKUgP3GLhjra4Hp2U53hsBlk4Hu0zHuKOytRf74BiX+/GZphfjVrArpb6nlDqplDqilHpWKeXf3edFQCuYXPIYiDssj7X2AbDMiyP2/OQYIsNM8gDVTW+esqIUbClJNbqUeXN3xP4yUKa1XgFUAV9zvyQhZjZ1cHMgLnmsbO0lOyGKZEuE1+5hNikWpsdSLZuULmpntZXlD7zEXY/tY2zi7ANK3qyyUpadQEpspEHVuc+tYNda/11rPfWk5h0g1/2ShJjZ1JLH+gCciqls6WWpF0frU0rSY2WO/SIcTs03njtG/8gEr53s5MXKtumP9Q6Pc7CpJ6CnYcCzc+x3AS+e74NKqbuVUhVKqQqr1erB24pQEohLHofGJqjrGmRZtvdXWJRkxNHaOxKQD5h95a9H26izDvLTj68hPzmG3+5pnP7Y4+804HBqrloSmMscp1w02JVSryilKmf4tf2M9/wzMAE8cb7raK0f0VqXa63L09IC+7uhMI6rfW9gBfuJtj60xqtLHaeUpLt6xtTKoRvn9fu9jeQlR3PNskxuW5/PntN2ajoHqO8a5EevVnPNskzWeGF3sC9dtKuN1vrqC31cKXUHcAOwVWutPVSXEDMqTrNgqxije3CMJC/OV3vSuw9OfTNiB9eSx1V5spbhvZrsQ7xda+O+baWYTIpby3P54StV/NPj+5lwaiLMJh64cZnRZbrN3VUx1wD3AzdqrQNv4lMEnMWTG0ZOtPcZXMnsHWvpIykmnKyEKK/fKz85hogwkzxAncHQ2AQPPHcMpeDmta7HgamxkXznQ8uxDY5Rbxvk2x8sI9MH/528zd0+lD8BIoGXJ4/5ekdr/Rm3qxLiPJZkuUakJ9v6uWRBYCxHO9Hex9LseI8dhXchZpNiQZo8QH0vp1N7jyQ+AAAUyElEQVTzmccPsKvaygMfWEZOYvT0xz60JpdtSzPoG5k46/VA5lawa60XeqoQIWYjLTaSFEsEJ9oCY8TucGpOtfdz+8YCn92zJD2WA43dPrtfIHjucCs7qqz8+/ZlfHJT4Tkfj4sKJy7KOz18jCA7T0VAUUqxJCuek+2BMSJtsA0yOuFkUabvWr8uTI+luXuYkXGHz+7p7369u55FGXHcvsF332CNJMEuAs6SrDhOdfQz4XBe/M0GOzX5DWiJh4/Cu5CpjVyNdnnsBTAy7uBoSy9XLA7c3i9zJcEuAs7izHjGJpwBMWo/2d6PUlAyeXSdL0yfNhVgy0K95WhLL+MOTXmBZw8Q92cS7CLgbClJxRJh5v6nj/j9dENVRz+FKRaiws0+u+dUsAdiszRv2N/get6wJj90ln9KsIuAkx4fxQ8+sorjbX38ZneD0eVcUFVH//RB076SEBNOUkx4wO3Q9ZY9dTaKUy0B3ftlriTYRUB6/7JMLi9N4yev1/jt9vmRcQf1tiFDzswsSLHIiB0YnXDwTp2dSwO4U+N8SLCLgPXFrQvpHR7n78c6jC5lRnXWQRxOPb0b1JcKU2JkxA4caOhheNzBpQsl2IUICGvyk8hKiOLVk/4Z7I12V7AWp1p8fu+CFAutPcOMTvj3Mwhv21VjxWxSbFyQYnQpPiXBLgKWUorV+YmcaPPP1TFTyw3zkmN8fu8F6bE4teunhlC2q7qLVXmJxAfR5qPZkGAXAa00I4562yDDY/43Mm20D5EQHU5CtO9DZdHk9E8on6bUPTjGkZbegD4Jab4k2EVAW5QRh9ZQ44dNr5rsw+QlG9N7pCjVQphJhXSwv11rQ+vAPuJuviTYRUArndyqf8oPA6zJPkS+AdMwABFhJgpTLVSF8MHWu2qsxEWGsTI3dNavT5FgFwGtYLJNrb+NTB1OTXP3sCHz61MWZcSFbJdHrTU7q7vYuCCFMHPoxVzofcYiqISZTSxMi53uyeIvOvpGGHM4DRuxg6uNQYN9yC+fP3jb6a5BmruHuSzAzy6dLwl2EfAWZcb53Yh9akWMkcFeOvn8odYaetMxb1a5zlW+LATn10GCXQSB0ow42npH6B32nx2oTX4S7IDf/TTjbRMOJ88faaMkPZaCFN/vIfAHEuwi4C3OfPecT3/RZB/CpCDbwBN5ClNiiDCbqOr0n6+Lt/3zs0dZ+M8vsr+hm4+syzO6HMNIsIuAN7Uyxp/a+Dbah8hKiCbcwAd3YWYTxWkWqkNkZczLxzt4Yk8jq/MTufuyYp+eWuVv3D3zVAjDZSdEERsZ5lfz7I0GLnU8U2lG3HTb2mD36K7TFKTE8NQ9mwz9huoPQvuzF0FBKUVphn+tjGm0D/tJsMfS0jPM4OiE0aV41dDYBBUNdt6/LDPkQx0k2EWQmFoZo7U2uhSGxiboGhglP8X4YJ/qLFnthztzPWnPaTvjDh2Su0xnIsEugsKijDi6h8bp6Bs1uhSau4cByE0y7sHplOmeMX7004w37KzqIiLMxLrC0Dn+7kIk2EVQKMtJAKCypdfgSqDRZvxSxyl5yTFE+uHOXE/bVWNlQ1GyT48g9GcS7CIoLMmKRymobPWDYPeDNexTzCbFwvRYqoJ4Kqa9d4SqjoGQO0zjQjwS7EqpryiltFJKvrLCEJbIMIpTLVS29BldCo32IWIizCRbIowuBYAFabHUBfHu0101XQAhd/zdhbgd7EqpPOB9QKP75Qgxf2U5CRzzgxH7VFdHpZTRpQBQnGahpWeYkfHg7Bnzdk0XKZYIlmTGG12K3/DEiP2HwP2A8csRREhbkZtIW+8IrT3DhtbRYB+iwA9WxEwpTotFa1djrGCjteatWlcXR5PJP76R+gO3gl0ptR1o0VofnsV771ZKVSilKqxWqzu3FWJGU0vddlQZ9+/L6dSG9mGfycK0WMA/DyNxV13XIB19o2xeINMwZ7posCulXlFKVc7wazvwdeDfZnMjrfUjWutyrXV5WlpottIU3lWSHktWQhRvnDIu2Nv7RhidcJLvR82nFqRbMJuUX23g8pQ9dXYANhbLMsczXbSlgNb66pleV0otB4qAw5NzibnAAaXUeq11u0erFGIWlFJctTidZw+2MDLuMGTp2+GmHgDKsv1nvjcyzMyCNAsn241/sOxp++rtpMZGUJTqP99I/cG8p2K01ke11ula60KtdSHQDKyRUBdGuqYsk6Exx3Q/bl/b39BNRJiJZdkJhtz/fBZnxnO8NbiC3eHUvFNno7wg2W8eVPsLWccugsrG4hQSosN58WibIfff39jNytwEIsL86/9aq/ISae0docUHD5a7Bkb56Rs1vH6ykweeO8a//qmS3bU2j9/n+SOttPWO8IGV2R6/dqDzWHfHyVG7EIYKN5vYtjSDlyrbGZ1wEBnmu+mYkXEHlS293HVpkc/uOVsbJueg95628cHVuRd9/+GmHr7/91MUpMTw9euWEBPhigqnU/Pj12r486EWJpyan358zfSuX4CxCSd3/HLvWfsJlILfvNPAr+9a77Gj6hyTdZRmxHJtWaZHrhlMpG2vCDo3rcrh6f3NfPqxCj66Po/rl2f55Ef1ypZexh2a8gL/e5C3ODOe+KgwdtdePNjHJpx86clDtPeOsKumi1dPdJIQHc7SrHj2nLbT0jPMmvxEGmxD3PbIO/zmHzawKi8RgB+/Vk1lSx//8cEy2npGuHxRGosz47j+R7v41vPH+csXLp3Xsw+nU7O/sZu4qDAWZ8bzwtE2ajoH+MnHVssyxxlIsIugs3lhCp+5fAFP72/m8789SMu1w9xz+QKv33eq7/ma/ESv32uuzCbF5oWp7KzuQmt9wW90zx9p5XTXIL+8Yx3hZhO/2FVHc/cwzxxsITM+iv/zweXctj6Ptt4RPvzwbj7x8z0kWsKJjwrnWGsft6zN5eMbzj7k4ps3LuPOx/bx/ZdO8S83LJ1T7cNjDu55fD87qqyYTYo/f24zP36tmpL0WK4ry5rX1yPYSbCLoKOU4qvXLub+9y/irl/t4yev13BreZ7Xt/hXNHRTlGohJTbSq/eZr8tL03ixsp2T7f0syTr/qp3XTnaSFhfJFYvSUEpxaUkqWmuOtfaxIC2W6AjXiDs7MZpHPlHOA88dIzYqjM7+EbaUpPLAjcvOueaVi9O5ZW0uj+9p4AtbS0iIDp9Vza09w3zlqcPsOW3ji1tLePjNWm748S4AfnybjNbPx7+e8AjhQSaT4mvXLmF03MlHH9mNtd97LX211hxo6GZNfpLX7uGurUsyMCku+GB5bMLJjiorl5emnTWqV0pRlpMwHepTlmbH89RnNvHoHet4/gtb+M2nNxAbOfN48Y5LChkZd/LsgeZZ1VtrHeDqH7xJRYOd792ykvu2lU7/5HVtWSbXL5fR+vlIsIugtigzjl/euY4m+zBfevKg1w7iaLANYRsco7zQf4M9LS6SzQtTeWJPI3/c38z2n+zirsf24XS++zXZVWOlb2SC65Z7/oFkWU4CZTnxPHuw5aLvHXc4+dozRwkzKV6573JuXut6LvClrSU89/nN/ORja2S0fgES7CLobV6Yyv3XLOKtGhsHGr1z/ufU/PraAv8NdoCvXruY/pEJvvKHwxxu7uW1k5386LVq/lbZxoTDyfNH2oiLCuPShd7ZHf6BFdkcbu6lwXbhvjWPv9PA3tN2vvGBZRScsYvXZFKsyE3ELKF+QRLsIiR8uDwPS4SZ3+9t8sr1KxpcKzam+rL4q2XZCfz585u5b1sp+//lai5ZkML/vFLNZx4/QNkDL/HMgRZuWJHltXX4169wTZ+8WHn+fYxaax57u571Rcl8aE2OV+oIdvLwVIQES2QYH1iZzZ8PtfKNG5eddx54vt6u7WJDUXJATA8syYqffnj64MfW8NzhVuyDY7T3jpCdGM2dlxZ67d65STGszE3gL4db+cx5ViodaOymwTbEF64qkR2l8yTBLkLGLWtz+f2+Jv5+rJ0Prbn4Jp3ZarAN0mAb4tN+uDHpYpIsEXzqkkKf3nP7qhz+/fnjVHX0Uzp5JuuZnjnQQlS4iWtk49G8yVSMCBlr8pPISYzmucOtHr3uzurJE3zkaLZZuXFVNmaT4v/fUXfOx5rsQzy9v5nrl2d7/KeqUCLBLkKGyaS4YWUWO6u76B0a99h1d1ZbyUmMlg6Ds5QaG8knNxXwh/3N/PKt09Ovf++lk9z04FsoBfe9r9TACgOfBLsIKe9bmoHDqdlR7ZnujxMOJ2/X2LisNFXmg+fg69ctYVVeIg+9WcuEw8n+BjsPvl5LdmI0T969iZzEaKNLDGgS7CKkrMpLIibCTEW93SPXO9zcQ//oBFtK5PCYuQg3m/jclQvp6BvltZOd/OyNOhJjwnnyno2szPO/lgyBRiaxREgxmxQrcxM5OHkghrvemTzBZ1NxikeuF0quXJRGRnwkX33mKPbBMe7dWjLdRVK4R0bsIuSszk/keGsfI+MOt691rLWX/OQYkrzchyYYhZlN/OOWYuyDY8RFhnHn5kKjSwoa8u1RhJzV+UlMODWVLb2UF7rXYreypY+yHP85Bi/QfOqSQqIjzGwoSiExRr45eoqM2EXIWT3ZVreiwb32Ar3D4zTah/zuGLxAEm428fENBSxM9+8du4FGgl2EnNTYSIrTLOw77d4D1KkzRM88QUgIfyDBLkLShqJk9tbbcTjn3+3xcLPrAexyCXbhZyTYRUhaX5RM/8gEp9r7532Ng43dFKTEeP0ADyHmSoJdhKT1Ra7liXtP2+b197XWHGzsYbWsuRZ+SIJdhKScxGhyEqPZO8+NSm29I3T2j7Laj09MEqFLgl2ErA1Fyew9bZ/XqUqHJjc4rZIRu/BDEuwiZK0rSqZrYIx629Cc/+7Bxm4iwkwXPBRaCKO4HexKqS8opU4qpY4ppb7riaKE8IUVua7VLEdbeuf8dw819VCWHe+1k4aEcIdb/yqVUlcC24GVWutlwPc9UpUQPlCaEUdEmImjzXPrGzPucHKkuVfm14Xfcne48U/Ad7TWowBa6073SxLCN8LNrqmUuY7YT7b1MzrhlPl14bfcDfZSYItSao9S6k2l1LrzvVEpdbdSqkIpVWG1eqYXthDuWp4Tz7GWPpxz2Kh0qMnVimCqNYEQ/uaiwa6UekUpVTnDr+24moglAxuB/xd4Sp3ntAGt9SNa63KtdXlamvSuFv5heU4C/aMTNNhn/wB1d52NzPgoOQxC+K2LdnfUWl99vo8ppf4JeEa71ovtVUo5gVRAhuQiICzPcY26jzT3zOpou5FxBzururhueZacmCT8lrtTMX8CrgRQSpUCEUCXu0UJ4SuLMuOIiwpjd+3sdqC+cLSN/tEJblqd4+XKhJg/d/uxPwo8qpSqBMaAT+n57PYQwiBmk2JjcQpv1c5uPPJURRMFKTFsLHavj7sQ3uTWiF1rPaa1vl1rXaa1XqO1fs1ThQnhK5sXpNBkH6Zxho1KddYBeofHAWiwDfJOnZ1b1uTKNIzwa3KCkgh5ly9Kh78c5+EdtfzbB5YSGWYGoKLezq0P7yY2MozNC1Jp6x0mzKT48Lo8gysW4sIk2EXIK0q18NF1eTyxp5GWnmEeu3M9Wmu+/dcTmJVicWYce+vtTDicfPumMjLio4wuWYgLkmAXAvjPDy0nKyGaH75SxTt1NmwDYxxq6uG7N6+YHqFrrWUKRgQECXYhAKUU91xezON7GnjguWMMjTlYlBHHzWtzz3qPEIFAOhgJMSkq3MyXri7hZHs/jfYhvrl9GWaThLkIPDJiF+IMH1ufT2ffKCtyE9hYnGJ0OULMiwS7EGdQSvHlbaVGlyGEW2QqRgghgowEuxBCBBkJdiGECDIS7EIIEWQk2IUQIshIsAshRJCRYBdCiCAjwS6EEEFGGXEuhlLKCjTM86+nEnqnNMnnHBrkcw4N7nzOBVrrix4abUiwu0MpVaG1Lje6Dl+Szzk0yOccGnzxOctUjBBCBBkJdiGECDKBGOyPGF2AAeRzDg3yOYcGr3/OATfHLoQQ4sICccQuhBDiAiTYhRAiyARUsCulrlFKnVJK1Silvmp0Pd6mlMpTSr2ulDqulDqmlLrX6Jp8QSllVkodVEo9b3QtvqCUSlRKPa2UOqmUOqGU2mR0Td6mlPry5L/pSqXU75RSUUbX5GlKqUeVUp1KqcozXktWSr2slKqe/N8kb9w7YIJdKWUGHgSuBZYCtymllhpblddNAF/RWi8FNgKfC4HPGeBe4ITRRfjQ/wJ/01ovBlYS5J+7UioH+CJQrrUuA8zAR42tyiseA655z2tfBV7VWpcAr07+2eMCJtiB9UCN1rpOaz0G/B7YbnBNXqW1btNaH5j8fT+u/8PnGFuVdymlcoHrgZ8bXYsvKKUSgMuAXwBorce01j3GVuUTYUC0UioMiAFaDa7H47TWOwD7e17eDvxq8ve/Am7yxr0DKdhzgKYz/txMkIfcmZRShcBqYI+xlXjd/wD3A06jC/GRIsAK/HJy+unnSimL0UV5k9a6Bfg+0Ai0Ab1a678bW5XPZGit2yZ/3w5keOMmgRTsIUspFQv8EfiS1rrP6Hq8RSl1A9Cptd5vdC0+FAasAX6mtV4NDOKlH8/9xeS88nZc39SyAYtS6nZjq/I97Vpr7pX15oEU7C1A3hl/zp18LagppcJxhfoTWutnjK7HyzYDNyql6nFNtV2llHrc2JK8rhlo1lpP/ST2NK6gD2ZXA6e11lat9TjwDHCJwTX5SodSKgtg8n87vXGTQAr2fUCJUqpIKRWB62HLcwbX5FVKKYVr7vWE1voHRtfjbVrrr2mtc7XWhbj++76mtQ7qkZzWuh1oUkotmnxpK3DcwJJ8oRHYqJSKmfw3vpUgf2B8hueAT03+/lPAn71xkzBvXNQbtNYTSqnPAy/heor+qNb6mMFledtm4BPAUaXUocnXvq61fsHAmoTnfQF4YnLAUgfcaXA9XqW13qOUeho4gGvl10GCsLWAUup3wBVAqlKqGfgG8B3gKaXUp3G1Lv+wV+4tLQWEECK4BNJUjBBCiFmQYBdCiCAjwS6EEEFGgl0IIYKMBLsQQgQZCXYhhAgyEuxCCBFk/i9Wy9kfpDsNDQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "N_sim = 10000\n",
    "e = np.random.randn(N_sim)\n",
    "te = np.arange(N_sim) * Ts\n",
    "\n",
    "_, d, _ = control.forced_response(Hu, te,e)\n",
    "plt.plot(te,d)\n",
    "np.std(d)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
